Global Attractors of Sixth Order PDEs Describing the Faceting of Growing Surfaces

A spatially two-dimensional sixth order PDE describing the evolution of a growing crystalline surface h(x, y, t) that undergoes faceting is considered with periodic boundary conditions, as well as its reduced one-dimensional version. These equations are expressed in terms of the slopes u1=hx\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_1=h_{x}$$\end{document} and u2=hy\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_2=h_y$$\end{document} to establish the existence of global, connected attractors for both equations. Since unique solutions are guaranteed for initial conditions in H˙per2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{H}^2_{per}$$\end{document}, we consider the solution operator S(t):H˙per2→H˙per2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$S(t): \dot{H}^2_{per} \rightarrow \dot{H}^2_{per}$$\end{document}, to gain our results. We prove the necessary continuity, dissipation and compactness properties.


Introduction
We study the long time dynamics of the following system arising in modeling of the self assembly of nanostructures, (1) Here, h : × R + → R is the height of a growing crystalline surface undergoing faceting during growth, δ > 0 is related to the deposition strength, W : R d → R is the anisotropy function derived from the surface energy anisotropy and ⊂ R d is the spatial domain with d = 1 or d = 2. Depending on the dimension we consider two specific cases for W given below, see (2), (3). We restrict our attention to special geometries and work with = (0, L) d , d = 1, 2 and furthermore we assume periodic boundary conditions for h. We show existence of global attractors for the above system and we shall also explain, why this result seems the best we can hope for. The PDE was introduced by Savina et al., see [20]. It has been derived by invoking Mullins' surface diffusion formula [15], a normally impinging flux of adatoms to the surface and a strongly anisotropic surface energy formula. The reduced evolution equation is obtained by carrying out a long-wave approximation. The choice of periodic boundary conditions is realistic as the patterns of the nanostructures statistically repeat throughout the domain, which is much larger than the length-scales of interest. Numerical simulations imposing these kinds of boundary conditions show good agreement with the experimentally observed behavior of crystalline materials undergoing faceting and coarsening [7,20]. We also notice that the analysis on periodic domains is easy to transfer for the numerical analysis of simulation schemes based on trigonometric interpolation. Such collocation methods are applied frequently to such problems.
We consider the following anisotropy functions W . If d = 1, then we take yielding a double well potential. In the two dimensional case, a naive generalization of (2), i.e. W (F 1 , F 2 ) = 1 4 (F 2 1 + F 2 2 −1) 2 is not appropriate, if we want to model growing pyramids, see [20]. For this reason we deal with where α, β > 0 are anisotropy coefficients. Formula (3) gives a quadruple well that is responsible for the faceting of the growing surface in shape of pyramids with four preferred orientations and hence preferred slopes. A constant A may be chosen such that W is always nonnegative.
In two related works, [10, Theorem 1.1] and [11, Theorem 2.1], we proved the existence of global in time weak solutions to (1) with periodic boundary data. There were no size restrictions on the data.
In [10,11] we showed only exponential bounds on the growth of solutions which is not particularly suitable for studying long time behavior. We will find the remedy here and we will show existence of a global attractor of (1) for d = 1, 2. The destabilizing term does not give us much hope to establish convergence to an equilibrium state. However, if we had a Liapunov functional, then we could hope to use methods based on Łojasiewicz inequality to show convergence of solutions to a steady state, see [19].
Our plan is to study first the one-dimensional problem, so that we can develop ideas that are used later also in the more complex case. It turns out that the trick applied in [11,20] works very nicely. Namely, after differentiating (1) with respect to x we obtain a slope equation for the new unknown quantity u = h x , see (4). One advantage is that we obtain a new conserved quantity, L 0 u dx = 0. This will imply that the semigroup generated by 3 has an exponential decay. Another advantage is, the resulting equation is similar to the convective Cahn-Hilliard equation, which has already been analyzed to some extent. Equation (1) may be interpreted as a convective Cahn-Hilliard (CCH) type equation of higher order, hence we call it the HCCH equation. Note that it is the gradient system perturbed by a destabilizing Kardar-Parisi-Zhang type term |∇h| 2 .
Here, we use ideas from the theory of infinite dimensional dynamical systems [4,18] combined with the available results on convective Cahn-Hilliard equation, e.g. [1,3,12]. Eden and Kalantarov noticed, see [1], that the structure of the lower order convective Cahn-Hilliard equation permits to deduce bounds implying the existence of an absorbing set. The same method can be applied here. We deduce from it the existence of an absorbing set in the H 1 topology and we extend this result to H 2 . Showing its compactness in H 2 requires further improvement of the regularity of weak solutions. Once we have achieved this goal, we may conclude the existence of a global attractor, see [14,Theorem 1].
We notice that, if we take the gradient of (1) with respect to the spatial variables in the two dimensional case, then the resulting system, see (8), has the structure which permits to carry the calculations we did for the one-dimensional problem. Thus, we establish the existence of the global attractor for the corresponding system, which is the result of the gradient of (1), and we call it the slope system, u = ∇h. Finally, we deduce from this existence result, the existence of a global attractor of the original Eq. (1), see Theorems 4, 5, 6.
We proceed as follows. In the next section we recall the notion of weak solutions and the necessary facts from [10,11]. In addition we state the main results. In Sect. 3 we prove the existence of absorbing balls in H 1 for the one-dimensional problem (4). This is done with the help of ideas taken from [1]. We also show in this section the necessary auxiliary facts. In Sect. 4, we study the system, which is obtained by taking the gradient of (1) and we call it the slope system. Its advantage is that we can use exactly the same method, as in the one-dimensional case to show the existence of an absorbing ball in H 1 . Next section is devoted to the proof of higher order regularity and compactness in H 2 of the absorbing balls, we use the parameter variation formula for this purpose. This is done in both case d = 1 and d = 2.
Finally, we discuss the results and future plans in Sect. 6.

Properties of Solutions and Main Statements
Im fact, we treat in [10,11] existence of solutions in cases d = 1 and d = 2 differently. For the one-dimensional problem we switch to a new variable, the slope u = h x , i.e. we differentiate (1) with respect to x. The resulting problem is where In [11] we adopted a natural definition of a weak solution of the one-dimensional problem (4). In order to express it we introduce the following notation.
is a weak solution to (4) provided that it fulfills, (see [11]), with u(x, 0) = u 0 (x), where g(u) = 1 2 u 2 and T = × (0, T ). In fact, the first integral denotes the pairing between u t and the test function ϕ.
We showed the existence of such solutions: For the two-dimensional problem (1) and (3) Furthermore, the analysis of the two-dimensional problem gets simplified after we transform Eq. (1) to a system for the slopes, u = (u 1 , u 2 ) = (h x , h y ) and we study The derivation of this equation and specification of W will be carried out in Sect. 5. For the purpose of proving our main results we should recast Eqs. (4) and (8) in the terms of the semigroup theory. Continuity of S(t), t > 0 follows from results in [10,11]. The uniqueness theorems imply that the family {S(t)} t≥0 has the semigroup property. It remains to establish strong continuity of the family {S(t)} t≥0 in the one-dimensional case. This will be done in Sect. 3.2. On the other hand, strong continuity of S in the two-dimensional case of (1) has been already established in [10].

Proposition 3 Let us denote by S(t)u 0 the unique solution u(t) to
The use of the language of the semigroup theory does not imply that we need to re-prove our existence results exploiting the analytical semigroup theory, see [5] or [2]. If we tried this, then we would repeat estimates specific for these problems presented in [10,11]. However, we will need additional regularity estimates, which we will establish with the help of the constant variation formula see Sects. 2.3, 5.1.
Here are our main results.

generated by Eq. (8) with periodic boundary conditions has a compact global attractor.
Once we show these results we may address the question of the behaviour of the solutions to the original problem (1). We notice that one can easily recover a continuous function f from its derivative and its mean. Thus, the above results imply:

Theorem 6 The semigroup generated by Eq. (1) has a global attractor in H 3
per for d = 1 and d = 2.
We note that we have numerical evidence of the existence of such an attractor in the onedimensional setting. Figure 1 shows similar pictures of the evolution as in [11] for two values of δ. For large values of δ a strange attractor seems to be existent, see the time-space plot in (a) and one particular solution in (a'). For smaller values we numerically expect stationary solutions as in (b) and (b)', or traveling wave (time-periodic) solutions. Note that once the structures form, the solutions stay in anḢ 2 per ball as the theory predicts. One might hope to be able to prove that at least for small initial data the L ∞ norm of u stays roughly below 1, independently of the value for δ. Our analytical result, however, gives us information of different nature. We can take bigger initial conditions and still the absorption is in the  4) is shown together with the decrease of the norm of the discrete solution and the three phase spaces (u, u x ), (u, u x x ) and (u x , u x x ). However, these runs calculated with a pseudospectral method discussed elsewhere [11], are particular examples with special initial conditions, fixed domain length and deposition parameter. The theory establishes a general result. We would like to remark, that more simulation results, also for the two-dimensional setting, can be found in [20].

Tools of Dynamical Systems
We will use the methods of the infinite dimensional dynamical systems, see the books by Hale, [4], Temam, [22] or Robinson, [18]. However, we will use the theorem guaranteeing existence of a compact global attractor as stated in [14]. The general theory stipulates that S(t) : H → H is a semigroup, where H is a Hilbert space. Following [14], we recall the necessary notions.
Let us suppose C 1 , Here, we note that any absorbing set attracts bounded sets.
For a bounded set B ⊂ H , we define its ω-limit set by the formula A compact global attractor for S(·) is a maximal compact invariant set.
Theorem 7 (see [14,Theorem 1]) Let us suppose that S(·) has a compact attracting set K . Then there is a compact global attractor for S(·) and A = ω(K ).
Alternatively, we could establish first that S(t) : Our line of argument, however, is based on Theorem 7. This result will imply our Theorems 4 and 5 once we show its assumptions are fulfilled. For this purpose we need Proposition 3. We also have to show the existence of a compact attracting set K . This will be achieved in two steps. First, we will establish existence of an absorbing set in H 1 . Next, by application of a different method, the existence of an absorbing set in H 2 and its compactness will be proved.
Note thatḢ 2 per is the correct choice of space for the slope systems (4) and (8), because we could not work with solution operators acting on lower order spaces due to the lack of uniqueness.

The Integral Representation of Solutions
The energy estimates become more tedious in two dimensions. Therefore, we choose a different approach to prove the higher oder absorption.
Let us consider the 2-d system, (28), for u = (u 1 , u 2 ). We note that we have the formula that has been used for the surface equation for h in [10]. Here the exponential operator is defined by where the right hand side is the inverse Fourier transform, whilef denotes the Fourier coefficients. For more details we refer to the cited work. It turns out that we can derive the same formula for solutions of the one dimensional problem (4). Indeed, since we have a unique weak solution, we may apply the Fourier transform to both sides of (4). The knowledge of the parameter variation formula for ODE's yields, change from estimate to estimate, but does not depend on the initial condition. This quantity may rely on the domain length and the deposition related parameter, L and δ, respectively. Numbers whose actual value is needed for balances with other estimates are denoted by C j , where j is an integer index, and these numbers are fixed.
In the second part of this section we will show that the semigroup S(t) :Ḣ 2 per →Ḣ 2 per is indeed strongly continuous. This will be done by a series of a priori estimates of Galerkin approximations and passing to the limit.

Absorbing Ball in H 1
Consider the HCCH Eq. (4) with periodic boundary conditions on a domain = (0, L) and initial condition u(x, 0) = u 0 (x). We extend the analysis from [11] by showing that the solutions are in fact absorbed into a ball whose radius does not depend on the initial value's norm. To prove this result we will need to combine several estimates that we want to formulate as separate statements. Subsequently, we write the L 2 -norm as · = · L 2 (0,L) and the L 2 scalar product by (·, ·). Other norms are equipped with a corresponding subscript.

Lemma 8 Weak solutions to Eq. (4) with u
Proof Application of the integral operator (− ) −2 :L 2 →Ḣ 4 per to both sides of Eq. (4) yields with g(u) = u 2 /2. The regularity guaranteed by (6) implies that the expression δ Hence, u t belongs to the same space and it may be paired with the left-hand-side of (12). Next, integration by parts, rearranging and using the Cauchy inequality with yield, d dt , with an L dependent constant, we get d dt Choosing 0 = 1/δ, we obtain d dt with a constant C 1 = C 1 (L , δ). In fact, by noting that we obtain the estimate (11).
The result shown above can be used for proving the existence of absorbing sets inḢ 1 , therefore one needs to take care of the right hand side in (11). This is done in the following lemma.

Lemma 9 The inequality
is true for all weak solutions of (4) with u 0 ∈Ḣ 2 per .
Proof We use u as a test function in the same transformed Eq. (12), . This is again allowed due to the regularity property (6). Putting the convective term on the right hand side, we estimate it as Furthermore, by using u 2 ≤ 1 4 u 4 L 4 + 1 = 1 4 u 4 L 4 + C 4 we finally derive (15) with C 2 = C 3 + C 4 , where both C 3 and C 4 depend upon L. Now we are able to prove the existence of the first absorbing set. Proof We define the 'energy' Then, by adding 4C 1 times estimate (15) to (11), we obtain Here we added and subtracted a small fraction of E 1 ( > 0). A rearrangement yields The H 1 term does not make any trouble, as can be chosen arbitrarily small. Furthermore, as we can estimate the right hand side and put the terms back on the left hand side to balance them with the L 4 term. In this way we derive Choosing sufficiently small yields and Gronwall Lemma leads to Remark We now know that u 2 ≤ C, u x 2 ≤ C and u 4 L 4 ≤ C, for a constant C independent of the initial condition that is undershot after a transient time. Since this case is one-dimensional this result leads to a uniform L ∞ bound on u. Furthermore, it was neither necessary to impose any restrictions to the deposition related parameter δ nor to the domain length L to achieve the result.
By the same method we can establish the existence of an absorbing set in the H 2 topology, but the argument is more involved. Possibly, we may show its compactness. However, this is of no use in the two dimensional case. This is why we will use a more general tool capable of handling both dimensional cases simultaneously. However, the starting point is the specific estimate like (18).

Strong Continuity of S(·)
We need to show that in the one-dimensional case Eq. (4) generates a strongly continuous semigroup. Since the original argument in [11] is based on the Galerkin method applied to Eq. (4), we will use it here.

Proposition 11
Let us suppose that u 0 ∈Ḣ 2 per , then S(t)u 0 ≡ u(t) converges to u 0 in thė H 2 per topology, as t → 0 + , where S(t) is the semigroup operator defined by (4).
We need the Gagliardo-Nirenberg inequality that holds on bounded domains ⊂ R n , for n ≤ 3. It states where j/m ≤ a < 1 and 1/ p = j/n + a(1/r − m/n) + (1 − a)/q , and where c 1 and c 2 are positive constants. For p = ∞ the fraction 1/ p is interpreted as 0.
Lemma 12 Let us suppose that u 0 ∈ B ⊂Ḣ 2 per , where B is a bounded subset ofḢ 1 per . Then, weak solutions to Eq.
Proof We note that and we estimate each of the three terms separately. Using Gagliardo-Nirenberg inequality (19) with n = 1, j = 1, p = 6, m = 4, a = 1/3, r = 2 and q = 2 we deduce Now, the inequality a 2 b 2 ≤ a 6 /3 + 2b 3 /3 for positive a and b implies The first term can be estimated as before, for the latter we again apply Gagliardo-Nirenberg inequality (19). We set j = 2, p = q = 3, m = 4, r = 2, n = 1, a = 12/23 so that Finally, we use ab ≤ a p / p + b q /q for conjugate numbers p = 23/18 and q = 23/5. This yields the overall estimate The last term in (21) is just bounded by C u x x x 2 , so that we derived (20).
Proof It is sufficient to show that u t ∈ L 2 (0, T ; (H 1 per ) * ). If we know this, then Proposition 1 (b), Lemma 12 and Eq. (4) imply that We act at the level of Galerkin approximation u N , see [11]. We apply the integral operator (− ) −1 :L 2 →Ḣ 2 per to both sides of Eq. (4) to derive We test this equation by u N t and estimate Collecting the (− ) −1/2 u N t 2 terms on the left hand side yields Applying the estimate (20) and integrating with respect to t we get, Since the right hand side is bounded uniformly in N due to the existence result established for the HCCH equation, we can pass to the limit and conclude that indeed our claim holds.

The Slope System in the Two Dimensional Setting
For the purpose of analysis of the two-dimensional spatial domain, we rewrite Eq. (1) as a system of slope equations. The surface with height h over the reference plane depends on the domain size, it grows due to coarsening that leads to an increase of the average size of the evolving structures. The slopes have a more dissipative character as the anisotropy of the surface energy forces the slopes to stay at a certain level that is independent of the domain size. We write u 1 = h x , u 2 = h y and note that the function in (3) is now used with u 1 and u 2 as arguments. In the evolution equation we need to calculate the gradient of W with respect to its arguments u 1 and u 2 , we denote it here by D u , − u 2 and further note that div D u W yields the second order linear and nonlinear terms. The fourth order term in the same potential stems from a corner regularization in the extended surface energyW = W + ( h) 2 /2. Now we transform Eq. (1) to a slope equation, using the same notation as introduced above. For this purpose we take the gradient of both sides of (1). If we set u = u 1 u 2 = h x h y = ∇h, then we will arrive at It is obvious from (28) that any solution to this equation is a gradient. The advantage of (28) is that it has the same structure as in the one-dimensional setting. Due to this fact and the prescription of periodic boundary conditions we have and hence again work with Sobolev spaces with zero mean. We establish the first absorption property:

Lemma 14
There is a constant C > 0 such that, for any bounded set B ∈ (Ḣ 2 per ) 2 , there exists a timet B , such that, for any u 0 ∈ B and any t ≥t B , the following uniform bounds hold We consider now (− ) −2 to be the inverse operator of the bi-Laplacian 2 :Ḣ 4 per ⊂ L 2 →L 2 and apply it to our new transformed system, The last equality is based on the observation that for any vector field X ∈Ḣ 1 This becomes obvious, after application of the Fourier transform to both sides, As before, we can test (30) by u t . However, this time we integrate over a two-dimensional domain, and as we deal with a system, we add the two components together, where we write shortly (− ) −1 u t 2 = (− ) −1 (u 1 ) t 2 + (− ) −1 (u 2 ) t 2 and keep this notation for all norms with arguments that are two-dimensional vectors, i.e. u L p = (u 1 ) 2 + (u 2 ) 2 L p . We arrive at A series of integration by parts based on u t = ∇h t yields, Using the identity in (32) we derive d dt corresponding to (11) in the one-dimensional setting.
Using u as a test function in (30) and adding the components yields Similarly, as before we use the identity which can be derived by the same argument as above.
We estimate the two terms with the wrong sign as where C depends on the domain parameter L as we applied Young's inequality to u 2 · 1.
Overall we have derived Multiplying (33) by α/(9D 1 ) and adding to the above estimate yields d dt We define the energy and proceed analogously as in the proof of Theorem 10. Hence, once again Gronwall Lemma yields the existence of absorbing sets in H 1 . and hence (36) holds. We notice that the dimensionality of the problem does not intervene here.
We may now establish new results based on (36).

Lemma 16
There is a constant C > 0 such that, for any bounded set B ⊂Ḣ 2 per , there exists a time t B =t B + 1 > 0, such that, for any u 0 ∈ B and any t ≥ t B we have Proof It is sufficient to show the following bound for any α > 0 We note that indeed and proceeding as in the proof of (36), we conclude that Hence, it follows from (10) and (36) for s = 1 + α, where α > 0, where we pick t 0 =t B and D u W (u) is the cubic nonlinearity. We notice that Furthermore, we can estimate the other nonlinear term by The last estimate is a consequence of Lemma 14. The uniformity of the constants of the last two estimates also comes from the uniform absorption of bounded sets of the energy (35). Our claim follows for t ≥ t B =t B + 1.

Compactness of Absorbing Balls
Using Lemma 15 we do not only show the existence of absorbing sets in H 2 but also their compactness. Therefore we make the following key observation.

Proposition 17
There exist α > 0 and a constant C(α) > 0, such that for any bounded set B ⊂ (Ḣ 2 per ) d , d = 1, 2, there exists a time t B = t B + 1, such that, for any u 0 ∈ B and any t ≥ t B we have Proof We use again formula (36), this time we take p = 2 + α, we use (10) too. We get We also observe that due to Lemma 16 we have, Combining these estimates with (40) we conclude that our claim holds for t ≥ t B = t B + 1.
We may complete the proofs of Theorems 4 and 5 in one stroke. Proposition 17 yields compactness of an absorbing ball in H 2 topology. On the other hand we have already established the strong continuity of the semigroup S(t). Thus, an application of Theorem 7 finishes the proof. Now, we prove the final assertion. We transfer the above results to the problem expressed in terms of the shape h in (1).

Proof of Theorem 6
Exactly as in the slope system it is sufficient to show existence of a compact absorbing set in the H 3 topology. First, we notice that it is easy to reconstruct a function h : = (0, L) d → R, when it is given its derivative u = ∇h and the mean m = h. Indeed, we have the following formulas: in case d = 1, These two formulas and Theorems 4 and 5 imply existence of compact absorbing sets, hence existence of a global attractor in H 3 .

Conclusions and Outlook
We have established the existence of global attractors inḢ 2 per for the slope Eqs. (4) and (8). This enable us to show the existence of global attractors in H 3 per for (1) in the 1+1D and 1+2D settings. On the way, we showed that solutions to (4) and (8) enjoy further regularity. For the one-dimensional case we succeed in deriving proper uniform estimates by repeated application of Gronwall inequality. As we needed uniform constants for the estimates, the work may seem somewhat tedious at certain points, e.g. during the application of Gagliardo-Nirenberg's inequality. Because of its repeated application this approach is not feasible the two-dimensional setting. Instead we reconsidered the constant variation formula from our previous work [10] to improve the regularity result. It turns out, once this approach is understood, the semigroup ansatz seems more elegant for this problem.
We are content with the results obtained for the presented equations. They coincide with the observations made with the help of a pseudospectral numerical method in the previous work [11], though we were not yet able to show or negate the existence of stationary or traveling wave solutions, which have been discussed in this publication. As we are not able to find a Lyapunov function, we were not in the position to use approaches based on the Łojasiewicz-Simon inequality (e.g. [13,21]).
We do not know much about the ω-limit set, but as Fig. 1 has already indicated, we expect to have time-periodic or stationary solutions for smaller values of δ and a strange attractor for increased values of the deposition rate dependent parameter. Note that once the structures form, the solutions in this figure stay in anḢ 2 per ball as predicted. The numerical simulations suggest that at least for small initial data the L ∞ norm of u stays roughly below 1, independently of the value for δ.