Superexponential droplet fractalization as a hierarchical formation of dissipative compactons Sergey Shklyaev,1,2 Arthur V. Straube,3 and Arkady Pikovsky2 1Department of Theoretical Physics, Perm State University, 15 Bukirev St., Perm 614990, Russia 2Department of Physics and Astronomy, University of Potsdam, Karl-Liebknecht-Str. 24/25, D-14476 Potsdam, Germany 3Department of Physics, Humboldt University of Berlin, Newtonstr. 15, D-12489 Berlin, Germany (Dated: January 27, 2010) We study the dynamics of a thin film over a substrate heated from below in a framework of a strongly nonlinear one-dimensional Cahn-Hillard equation. The evolution leads to a fractalization into smaller and smaller scales. We demonstrate that a primitive element in the appearing hierar- 0 chicalstructureisadissipativecompacton. Bothdirectsimulationsandtheanalysisofaself-similar 1 solution show that the compactons appear at superexponentially decreasing scales, what means 0 vanishingdimension of thefractal. 2 PACSnumbers: 47.53.+n,68.15.+e,47.20.Dr,05.45.Df n a J Introduction.—A vast number of intriguing pattern- Ref. [8]), where h is the local thickness and 7 formation phenomena can be described with high-order 2 BMh2 nonlineardiffusionequationsofCahn-Hillardtype. Since f = Boh3+ , g =h3. (2) ] theirintroduction[1],theseequationshavebeensuccess- − 2(1+Bh)2 S fully applied to a great variety of natural and techno- P DimensionlessBond(Bo),Biot(B),andMarangoni(M) logicalprocessessuchasphase separationinbinary mix- . numbers determine the levels of the gravity, of the heat n tures, alloys, glasses, and polymer solutions (see, e.g., flux at the free surface, and of the convective flow, re- i surveys [2]), topology transitions in a Hele-Shaw cell [3], l spectively. Although function f here has a rather com- n dynamics of layered systems [4], thin films [5], competi- [ plex form, for h 0 one can set f 0.5BMh2. This tion and exclusion of biological groups [6], and aggrega- → ≈ approximation holds also for moderate values of h, pro- 1 tion of aphids on leaves[7]. In the thin film context, nu- vided the gravity can be neglected, Bo 0, the heat v mericalstudiesofanamplitudeequationofCahn-Hillard → 6 transfer at the free surface is poor (B small) while the type[8,9]haveevidencedfilmruptureleadingtothefor- 9 thermocapillary effect is strong (M large). The function mation of a cascade of “drops” and “fractal-like finger- 9 g(h), which can be referred to as “mobility”, is conven- 4 ing” [10] comprising the gaps or “dry spots” [9] between tionally non-negative, g(h) 0, as this prevents against . the drops. These findings have been supported by direct ≥ 1 the fast growth of the short-wave perturbations. 0 simulations of the Navier-Stokes equations [12]. Assuming the limiting form of f, after an appropriate 0 The goal of this paper is to describe this cascade as a rescaling of the time, the field, and the coordinate we 1 hierarchical formation of dissipative compactons. Com- arrive at our basic equation : v pactonisawell-knowncompact(i.e. withfinitesupport) i traveling-wave solution, which emerges in conservative h + h2h +h3h =0. (3) X t x xxx x systems with nonlinear dispersion [11]. Its stationary r Noteworthy, Eq. (3)(cid:0)is invariant un(cid:1)der the scaling a analogue with compact support appears in dissipative systems with nonlinear dissipation and, therefore, can h p2h, x px, t p−2t, (4) be referred to as a stationary “Dissipative Compacton” → → → (DC). Below we demonstrate that a DC is a primitive meaning that a thinner film evolves slower. element mediating the formation of hierarchical fractal Steady state and its stability.—We start our analysis structure, and characterize the fractal properties of this by considering positive stationary solutions h = H(x) structure quantitatively. of Eq. (3). Looking for symmetric patterns, after one integration we obtain We focus on a one-dimensional Cahn-Hilliard equa- tion describing dissipative evolution of a conserved field ′′′ ′ HH +H =0 (5) h(x,t) with primes denoting d/dx. Equation (5) admits a com- h +[f(h)h +g(h)h ] =0. (1) pactsolutionH(x) in the formofa DC ora “touchdown t x xxx x steady state” [13], nonvanishing for x l only: | |≤ In the context of the dynamics of a thin film over a 1 substrate heated from below, this equation describes a x= √π erf ln H , =maxH(x), (6) surface-tension-driven convection (see, e.g., Eq. (4) in ± H r2 H! H x 2 dissipative compacton a nonzero volume and can be destabilizing, if nonlinear 1 corrections are retained. l ) x Thus, although a DC is stable with respect to per- PSfragrepla ements ~H( 0.5 H turbations of zero volume, the instability is possible if the volume of the DC is changed, as confirmed by our 0 -2 -1 0 1 2 numerical simulations of Eq. (3). We detect the tempo- x ral decay even for finite-amplitude perturbations of zero FIG. 1: (Color online). The shape of the base DC, H˜(x). volume. For the perturbations changing the volume but keepingconstantthelengthoftheDC,wefindabreakup of the DC with the emergence of a complex structure. where erf(z) = 2/π ze−t2dt. Solution (6) presents Evolutionary problem.—We now study the formation 0 a self-affine one-parameter family of DCs parametrized of a fractal, hierarchical structure of DCs, illustrated in p R by . For a thin film, the DC describes the stationary Fig. 2. We discretize Eq. (3) in the computation do- H profile of a drop with the amplitude and zero contact main 0 x d with periodic boundary conditions, with H ≤ ≤ angle. Owing to scaling (4), any DC can be expressed the number of nodes N = 1000 and apply the Newton- in terms of the base DC H˜(x) having =1, see Fig. 1. Kantorovich method [8]. We choose a distorted uniform H Thus, the profile of a DC and half its length l obey: profile h(x,t=0)=1+acos(2πx/d) with a=0.1 as an initial condition. Results of computations are presented H(x)= H˜ x/√ , l=√π . (7) in Figs. 2 and 3. There we also compare the numer- H (cid:16) H(cid:17) H ically obtained profile h(x) having local maxima h(mn), The property of self-affinity is a necessary prerequi- n = 1,2,... with the DC profiles with = h(n), DC(n) m site for the emergence of fractal structure described by below, what indicates that the initial staHte develops into Eq. (3), as we discuss below. As it follows from Eq. (7), a hierarchical structure of DCs of different amplitudes. DCs become narrower for smaller amplitudes – contrary This fact allows us to increasethe efficiency ofthe nu- to other examples of compactons where typically the merics significantly: because after their formation the width is amplitude independent [11]. In a more general DCs remain constant in their bulk, we exclude these do- situation, when for small h the functions f(h) and g(h) mains from numerical simulations and impose the corre- scale like f/g hγ, no solutions are possible for γ 2 sponding boundary conditions for the still evolving do- ∝ ≤− [13]. Atγ >0,thesolutionsaresoliton-likebecausetheir mains between the formed DCs. Thus, while proceeding support is no more compact. The case γ = 0 results in to smallerDCs, we canconsiderablyrefinethe meshand constantl. Compactsolutionssatisfyingtherequirement also increase the time step. Therefore, we not only re- of l 0 as 0 exist in the rangeof 2<γ <0 only. solve high-order DCs with the accuracy consistent with → H→ − Thus, the consideredabovecase γ = 1 is the only inte- that at previous stages, but also maintain the computa- − ger index possessing self-affine compactons. tional efficiency. This strategy provides reliable results To explore the stability of a DC, we introduce a small up to n=4. perturbation ξ(x)exp(λt) of H(x), where λ is the Theobservedstructurealongwiththepropertyofself- ∝ growth rate. By linearizing Eq. (3), we obtain affinity suggests that the formation of higher-order DCs λξ+ H3 ξ′′+H−1ξ ′ ′ =0. (8) nreemveinrisstcoepnstaonfdtthheeCdrayntsoprotssetb.etHweereen,DaCDs,Cfoprlmayasfararcotalel h (cid:0) (cid:1)i of a primitive element, mediating the fractalization. To Assuming ξ( l)= 0, we multiply Eq. (8) by ξ′′+H−1ξ characterize properties of this fractal quantitatively, we ± and integrate by parts to arrive at an integral relation plotinFig.4(a)thevariationofL ,thedistancebetween n l H′′ 2 l ξ ′ 2 λ ξ′ ξ dx= H3 ξ′′+ dx, Z−l(cid:18) − H′ (cid:19) −Z−l "(cid:18) H(cid:19)# (9) 100 which is closely related to the variational principle for 10-2 Eq.(1)[14]andthefactthattheLyapunovfunctionalhas hh 10-4 alocalminimumontheDC.AsH 0,boththeintegrals 101 in Eq. (9) are non-negative and ≥the perturbations are 10-6 103102 nongrowing, λ≤0. This result, howevPeSrf,radgoerespnlao tegmueanrt-s 1 2 3 4 5 105104 t antee against the instability, as there exist two modes of x neutral stability, λ=0, satisfying ξ( l)=0. One mode, ± FIG. 2: Fragment of the evolution of the field illustrating ξ1(0) =H′,reflectstranslationinvarianceandcannotgive hierarchical formation of droplets. Notice logarithmic scales rise to instability. Another mode, ξ(0) =H xH′/2,has of the timeand thefield. 2 − 3 3 Self-similar solution.—To shed light on the hierarchi- cal formation of DCs and to alternatively support the t=10 2 t=40 conclusions about the fractal dimension and the super- h t=400 exponential scaling, we construct self-similar solutions, 1 t=8000 which originate from the rescaling property, Eq. (4). By (a) seeking the solution of Eq. (3) in the form PSfragrepla ements 0 0 2 4 6 8 10 h=t−1G(η), η =x√t, (12) x 0.6 0.03 (b) (c) we arrive at an ordinary differential equation for G(η): 0.4 0.02 h h ηG′ 2G+2 G2G′+G3G′′′ ′ =0, (13) − 0.2 0.01 where primes stand for(cid:0)d/dη. Numer(cid:1)ical solutions of 0 0 Eq. (13) with various initial conditions all demonstrate 2 4 6 8 2.8 3.2 3.6 a qualitatively similar behavior of G(η), which displays x x aninfinitenumberofoscillationsofincreasingamplitude. FIG. 3: (Color online). Evolution of the film profile for d = Twonumericalsolutionscorrespondingtodifferentinitial 10. Panels (b) and (c) are zoomed in fragments of panel (a). conditions are shown in Fig. 5. LinesrepresentnumericalresultsaccordingtoEq.(3),circles For large G, where we can estimate d/dη ǫ1/2 with show the profiles of corresponding DCs as in Eq.(6). ǫ G−1 1, the first two terms in Eq. (∼13) become ∼ ≪ negligible in comparisonwith the last two terms. In this theneighboringDCsofn-thand(n 1)-thorders,versus limit, Eq. (13) is reduced to Eq. (5) with G(η) instead thebase2l ofDC(n). Thenumerica−lresultsfordifferent of H(x). Therefore, G(η) can be approximated by the n solutionfor a DC [seeEq.(6)andthe insetin Fig. 5] ev- d fit well a power law: erywhere except for its tails, where G is no longer large. L α(2l )β, α 0.2, β 1.25. (10) As a result, G(η) looks like a sequence of DCs with su- n n ≈ ≈ ≈ perexponentially growing amplitudes G exp[Ak] and k Note that deviations from this law for the points related ∼ widths ∆η = 2√πG ; the positions η of maxima for k k k to the biggest DCs stem from the initial condition. On large k satisfy η √πG , see markers in Fig. 5. k k theotherhand,forhigherorderstheself-similarityofthe ≈ To specify the superexponential growth of G with k, k formation of DCs is evident from Fig. 4(a). we construct a mapping G G valid for large G. k k+1 Because β > 1 in Eq. (10), with the increase in n the In the range of η η < √→πG , G G H˜(x ) with k k k k ratio L /l diminishes implying that the smaller daugh- | − | ≈ n n x (η η )/√G . To bridge the solution for DC(k) ter DCs tend to occupy the whole space between their k ≡ − k k with that for DC(k+1), we substitute a representation bigger parent DCs. The fraction of dry spots tends G = ε−2ζ(y), y = (η η √πG )ε, ε G−1/6 into to zero and, therefore, the fractal dimension of this − k − k ≡ k Eq. (13) and neglect the terms ε, what yields set equals zero. Furthermore, for large n we can ne- ∝ glect the distance between DC(n) and DC(n+1) and put y ζ′+2 ζ2ζ′+ζ3ζ′′′ ′ =0. (14) L 2l . As a result, Eq. (10) entails a remarkable 0 n n+1 ≈ superexponential scaling of ln with n: Here,primesdenoted/d(cid:0)yandy0 =η(cid:1)k/√Gk+√π 2√π. ≈ As Eq. (14) admits no analytical solution, we solved it log(l ) βnlog(l ). (11) n ∝ 0 numerically,withtheinitialconditionaty = ε−2(√π x): dpζ/dyp = ε2(p−2)dpH˜/dxp with p = 0−,1,2,3 an−d 1 20 dpH˜/dxp evaluated via Eq. (6) at = 1, which ensures (a) (b) H PSfragrepla ements logLn10 -- 021 d =8 logG k1011 268 wtthheeetmmakaaettccihhniitnnoggawwcicittohhutnthhteetdhgearctoawHy˜iinn′′gg=ttaa−iillloonffHDD˜CC−((kk1)+.[1cT)fo.atpEeyqr.≫fo(r5m1)], d=10 and obtain -3 4 -1.5 -1 -0.5 0 0.5 1 2 3 4 5 6 G =G1/3ζexp(ζ′′+1). k+1 k log102ln log10Gk+1 By determining ζ, we end up with the transformationof FIG. 4: (Color online). (a) The distance Ln between two neighboringDC(n) andDC(n−1) versusthebase2ln ofDC(n). Gk →Gk+1 [see Fig. 4(b)], which fits well a power law Squaresandcirclesarenumericalresultsford=8andd=10. G 40G2.83. (15) Dotted line is a fit, Eq. (10). (b) Mapping logGk+1(logGk) k+1 ≈ k calculated in the framework of Eq. (14) (circles); dotted line Equation(15)showsasuperexponentialgrowthforG k corresponds to the asymptotic law, Eq.(15). with k, which is expected as required by the self-affinity 4 10 amplitude equation and, alternatively, by constructing a 8 4⋅103 self-similar solution, we show that this structure of DCs G 2⋅103 is afractal,characterizedby superexponentiallydecreas- G 6 ing amplitudes and lengths of smaller droplets, and thus 10 0⋅100 g 4 100 225 350 having zero dimension. The dissipative compacton is a o l (cid:17) primitive element mediating the fractal structure com- PSfragrepla ements 2 prising the dry spots between the compactons. It should 0 be also noted, that a number of effects, such as inter- 0 1 2 3 4 molecular interaction between liquid and solid, contact log10(cid:17) angle dynamics, evaporation, etc. become of crucial im- FIG.5: (Coloronline). FunctionG(η)insolution(12). Initial portance, when the free fluid surface touches the solid. conditions at small η are G ≈ 0.5η2 (solid line) and G ≈ An extension of the theory above that includes these ef- 1−0.1η2 (dashed line). Squaresand triangles show thelocal fects remains a challenge. maximaGk ofG(η);forlargek,themaximaGk approachthe lawG=η2/π (dottedline). Inset: acomparison ofapieceof Acknowledgements.—We are grateful to A. Nepom- G(η) with a single DC (circles), Eq. (6). nyashchy, A. Oron, M. Zaks, Ph. Rosenau, D. Lyubi- mov, and D. Goldobin for stimulating discussions. The research was supported by German Science Foundation andthesimilarbehaviorforthelengths,seeEq.(11). As (projectsNo.436RUS113/977/0-1andNo.STR1021/1- the amplitude l2 [cf. Eq. (7)], the exponent 2.83 in 2) and Russian Foundation for Basic Research (project ∝ Eq. (15) is in reasonable agreement with 2β in Eq. (10), No. 08-01-91959). obtainedwithinthe evolutionaryproblem. The factthat the correspondence is not perfect is not surprising as Eq. (15) is the asymptote of extremely large t (i.e. large k), while Eq. (10) is a fit obtained for the early stage of [1] J. W. Cahn and J. E. Hilliard, J. Chem. Phys. 28, 258 the evolution (small k). Nevertheless, we see that the (1958). self-similar solution is closely related to the hierarchical [2] A. Novick-Cohen and L. A. Segel, structure of DCs described by the evolutionaryproblem. Physica (Amsterdam) 10D, 277 (1984); Mathemat- Finally, we stress that the relation between the self- ical methods and models in phase transitions, edited similar solution and the spatially periodic solution as in by A. Miranville (Nova Science Publishers, New York, Fig. 3 is not simple. The whole structure of successive 2005), Chap. 3. DC-like solutions, h(x,t), obtained via Eq. (13), moves [3] R. E. Goldstein, A. I. Pesci, and M. J. Shel- with the time toward the point x = 0, whereas DCs ley, Phys. Rev.Lett. 70, 3043 (1993), ibid. Phys. Fluids 10, 2701 (1998). H(x), which are born as a result of evolution according [4] D.V.LyubimovandS.V.Shklyaev, Fluid Dyn.39, 680 to Eq. (3), are stationary objects. However, the long- (2004). time evolutions of both these solutions show the sim- [5] A. Oron, S. H. Davis, and S. G. Bankoff, ilar displacement of the gaps between DCs by higher- Rev. Mod. Phys. 69, 931 (1997), T. G. Myers, order DCs. This argument becomes transparent, if we SIAMRev. 40, 441 (1998); R. V. Craster and observe the self-similar solution “stroboscopically”. Let O. K. Matar, Rev.Mod. Phys.81, 1131 (2009). us consider a self-similar solution at moments of time [6] D. S. Cohen and J. D. Murray, J. Math. Biol. 12, 237 (1981). t =η2/x2. Thecorrespondingfieldprofile(12)describes k k 0 [7] M. A.Lewis, Theor. Popul. Biol. 45, 277 (1994). the formationofDCs upto the k-thorderinthe domain [8] A. Oron, Phys. Fluids 12, 1633 (2000). 0 x<x0+∆ηk/√tkwithDC(k)centeredatx=x0. As [9] S.J.VanHooketal.,J. Fluid Mech. 345, 45(1997),ibid. ≤ thegrowthofGk withkissuperexponential,thehighest- Phys. Rev.Lett. 75, 4397 (1995). order DC dominates the pattern, what ensures that the [10] L. Y. Yeo, R. V. Craster, and O. K. Matar, fractal made of the dry spots has zero dimension. Phys. Rev.E 67, 056315 (2003). Conclusions.—We haveappliedtheconceptofdissipa- [11] P. Rosenau and J. M. Hyman, Phys. Rev.Lett. 70, 564 (1993); P.Rosenau, Phys. Rev.Lett.73, 1737 (1994). tive compactons to the evolution of a thin film within [12] W. Boos and A.Thess, Phys.Fluids 11, 1484 (1999). a framework of the generalized one-dimensional Cahn- [13] R. S. Laugesen and M. C. Pugh, Hilliard equation. We have shown that as a result of a Eur. J. Appl.Math. 11, 293 (2000). rupture, the thin film evolves into a hierarchical struc- [14] R.S. Laugesen and M. C. Pugh, J. Diff. Eqns.182, 377 ture of drops, which can be represented by dissipative (2002);L.Pismen, PatternsandInterfacesinDissipative compactons of different scales. By efficiently solving the Dynamics (Springer, Berlin, 2006).