Emergence of foams from the breakdown of the phase field crystal model
EEmergence of foams from the breakdown of the phase field crystal model
Nicholas Guttenberg , Nigel Goldenfeld and Jonathan Dantzig Department of Physics, 1110 West Green Street, Department of Mechanical science and Engineering, 1206 West Green Street,University of Illinois at Urbana-Champaign, Urbana, Illinois, 61801-3080.
The phase field crystal (PFC) model captures the elastic and topological properties of crystals witha single scalar field at small undercooling. At large undercooling, new foam-like behavior emerges.We characterize this foam phase of the PFC equation and propose a modified PFC equation thatmay be used for the simulation of foam dynamics. This minimal model reproduces von Neumann’srule for two-dimensional dry foams, and Lifshitz-Slyozov coarsening for wet foams. We also measurethe coordination number distribution and find that its second moment is larger than previously-reported experimental and theoretical studies of soap froths, a finding that we attribute to thewetness of the foam increasing with time.
PACS numbers: 47.57.Bc, 64.70.D-, 47.54.-r
Computational modeling of large-scale systems usu-ally involves either detailed molecular dynamics simula-tions, in which every particle must be tracked, or a highlycoarse-grained model in which the underlying symmetriesare known and are used to derive a continuum modelfrom the underlying physics. Molecular dynamics is lim-ited to small systems and/or to very short times, andcoarse-grained models tend to fail at places where theunderlying symmetries are broken, such as at defects anddislocations. The phase field crystal (PFC) model [1] isan intermediate approach with the advantage of diffusivetime scale, but with atomic scale resolution of molecu-lar dynamics. The PFC model can be used to capturethe dynamics of defects and dislocations in large crystals[2, 3], the dynamics of grain interactions [4], moleculardynamics of vacancies[5] and even nonlinear elasticity[6].In addition, it is amenable to methods such as coarse-graining and adaptive mesh refinement [7].The dimensionless phase field crystal equation [1] ∂ψ∂t = ∇ [( ∇ + 1) ψ + rψ + ψ ] (1)is a density functional theory [2], best thought of asarising from phenomenological and symmetry consider-ations. In particular, this is the simplest class of modelsappropriate for systems whose dynamics is governed byminimizing departures from periodicity[8], as opposed tothe situation in other materials processes, such as spin-odal decomposition, where the dynamics minimizes de-partures from spatial uniformity. In principle this ap-proach only holds for small values of the undercooling α ≡ − r , beyond which the strong nonlinearity may inprinciple overwhelm the crystal’s symmetries, leading tosuch artifacts as merger or dissolution of the ‘atoms’ ofthe model.In this Rapid Communication, we explore an interest-ing feature of the phase field crystal model in the limitof large undercooling, which we show turns out to facili-tate a simple, scalar and minimal model of foams. In the appropriate density regimes, the behavior of the phasefield crystal equation at large undercoolings is that theatoms of the crystal lattice begin to merge. However, theinterstitial spaces between the atoms are preserved andbecome line solitons. The result is a coarsening foam-likestructure. We analyze the process by which this insta-bility in the equation of motion occurs, and use this un-derstanding to propose a minimal, modified PFC modelthat is capable of describing quantitatively both wet anddry foams at the level of a continuum, scalar theory thatis computationally efficient. Equilibrium Phase Diagram:-
The phase diagramof the PFC model has been computed for smallundercoolings[1]. We shall use the same methods to con-struct the phase diagram at larger values of the under-cooling in order to see what may be found there. First,however, we will convert the PFC equation and PFC en-ergy to a nondimensionalized form with respect to theequilibrium liquid (constant phase) density. The energyminimizing density for the constant phase is ψ = ±√ α ,so we will introduce a nondimensional order parameter φ ≡ ψ/ √ α . This gives us the following free energy: F = (cid:90) V φ ( ∇ + 1) φ + α (cid:18) − φ + 14 φ (cid:19) dV (2)and corresponding equation of motion: ∂φ∂t = ∇ (( ∇ + 1) φ + α ( − φ + φ )) (3)We then construct the phase diagram at large α bycalculating the energy minima of the one-mode approxi-mations for the constant phases L ± : φ = φ , where thesubscript ± refers to the sign of the average density φ ;the striped phase S : φ = φ + A S cos( x ), where A S = (cid:112) α (4 − φ ) /
3; and the two triangular lattices ∆ ± : φ = φ +( A ∆ ± B ∆ )(cos( (cid:126)k · (cid:126)r )+cos( (cid:126)k · (cid:126)r )+cos( (cid:126)k · (cid:126)r )), wherethe (cid:126)k j are the lattice vectors of the regular triangular lat-tice, A ∆ = − αφ / B ∆ = (cid:112) α (15 − φ ) / a r X i v : . [ n li n . PS ] M a r Even though the one-mode approximation is not accu-rate at these large values of α , we use this approximationto understand heuristically the behavior observed, thusmotivating our form for the modified phase field crystalmodel introduced below. Our main results, for the mod-ified phase field crystal model, are fully time-dependentand independent of this approximation. FIG. 1: (Color online). Phase diagram associated with thefree energy in Eq. (2). Invariant reactions occur at α = 2 . α = 5 .
45. No further topological changes occur for α > The equilibrium phase diagram is obtained by substi-tuting the various ansatz into Eqn. (3) and performinga common tangent construction. The result is shown inFig 1. We find a series of previously undiscovered invari-ant points as α increases. At α = 2 .
56, the coexistenceregion between the striped and triangular phases disap-pears, giving rise to a coexistence between stripes andone of the liquid phases. Then, at α = 5 .
45 the stripedphase vanishes, leaving only a region of immiscible liquid-liquid coexistence, and the regions comprised entirely ofthe liquid phases.This can be understood by considering the influenceof the two terms in the free energy. At small α , thewavelength selection term that was added to constructthe PFC equation is dominant. As α increases, the localterm that drives φ to have the values of ± ∇ ) ) φ is a small perturbation and the equi-librium phases are determined by the average value of φ / − φ /
4, which is minimized by the constant phases.
Dynamics:-
The phase diagram we have just computedrepresents the behavior of this system in its final equilib-rium state. However, the approach toward that equilib-rium changes as α increases. Starting from a hexagonallattice, the system attempts to coarsen into two regions,one composed of L + , the other composed of L − . How-ever, because of the small residual wavelength selection,there is a finite energy cost to removing spatially pat-terned structures. The consequence of this is that the PFC ‘atoms’ are dynamically conserved by this infinites-imal energy barrier, which halts the coarsening dynamics.If the average density of the system is shifted, e.g., by de-liberately reducing φ over time, this can provide enoughenergy to destabilize the PFC atoms, but stripe-like cellwalls will still remain stable for a longer period as theycan move perpendicular to their length to create largebubbles of one of the liquids. FIG. 2: a. Final state obtained from Eqn. (3) with α = 16and φ = − .
41 after coarsening. Residual atoms coexist withlarge foam cell regions. b. Simulation using the modified PFCequation for α = 20, quenched with φ = 0 . φ = − .
6. The result is a dry foam structurewith no residual atoms. c. Simulation of the modified PFCequation for α = 20 and φ = 0 .
25, allowed to coarsen withoutdraining. The result is a wet foam structure with circularbubbles.
The end result of this is that a foam-like state appears(Fig 2a), which coarsens to a stationary point at whichthere are no more free adjustments that can be madeto approach the equilibrium state. Thus, even a smallperturbative addition of wavelength selection is sufficientto stabilize a foam-like structure.Although this behavior is interesting from the point ofview of understanding the properties of the phase fieldcrystal model, it does not behave like a physical foam.The arrested coarsening and inability to get rid of resid-ual atoms prevent this from being used as a model forstudying coarsening foams. We now alter the free energyso as to destabilize the residual atoms while retaining theoverall cellular structure, so that we may recover physicalfoam dynamics.In effect what we wish to do is to weaken the wave-length selection while encouraging k = 0 structures. Astripe has one k = 0 direction and one direction with theselected wavelength, whereas an atom has two directionswith the selected wavelength. If the selected wavelengthsits at a local energy minimum, but the global minimumis at k = 0, then this should help remove residual atomsin favor of bubble interfaces.The PFC wavelength selection operator ( ∇ + 1) cor-responds in k -space to E ( k ) = (1 − k ) (4)We modify this by introducing two parameters k and bE mod ( k ) = ( k/k ) (1 − ( k/k ) ) + b ( k/k ) (5)and then choose k so as to retain minima at k = ± k = [3 / (2 + √ − b )] / . The modified form of the freeenergy in real space becomes F mod = (cid:90) V φ (cid:32) k ∇ (cid:18) k ∇ + 1 (cid:19) − bk ∇ (cid:33) φ + α (cid:18) − φ + 14 φ (cid:19) dV (6)We can then vary b to control the relative depths of theminima at k = ±
1. The most extreme effect is achievedat b = − /
3, at which the minima at k = ± b throughout the rest of this paper.The dynamics of the modified PFC equation give riseto the structures shown in Figs. 2b and 2c. If the sys-tem is quenched from a disordered state, one recoversLifshitz-Slyozov coarsening dynamics [9] correspondingto the physical behavior of a very wet foam. On theother hand, if the system is drained from a state with pos-itive average density to one with negative average den-sity, polygonal cell walls form and coarsen by topologicalrearrangements. Results:-
We measured the coarsening dynamics andstatistics of the eventual scaling state of the modifiedPFC equation in order to compare with physical foams.In order to measure the coarsening we are interested inthe quantity (cid:104) r (cid:105) ( t ), the average bubble radius. We mayeasily count the total number of bubbles and so determinethe average bubble area as (cid:104) A (cid:105) ( t ) = A total /N ( t ). If weassume that deviations from circular geometry are small,then we may approximate (cid:104) r (cid:105) ( t ) ≈ (cid:112) (cid:104) A (cid:105) ( t ). FIG. 3: Coarsening behavior of wet and dry foams in themodified PFC equation with α = 20, compared with the cor-responding theoretical coarsening laws. For a wet foam, we expect the bubble growth to pro-ceed as diffusive grain growth, so that (cid:104) r (cid:105) ( t ) ∝ t / . Asan example, we simulate a 1024x1024 system with α = 20and quench into an average density φ = 0 .
25. We ob- serve a scaling (cid:104) r (cid:105) ∝ t . (Fig 3, filled circles), consistentwith the predicted growth law.For a dry foam in two dimensions, von Neumann’s lawfor bubble area growth [10] is dAdt = κ ( n −
6) (7)where κ is an effective diffusivity, and n is the coordina-tion number of the bubble. This implies that the averagebubble radius should scale as (cid:104) r (cid:105) ( t ) ∝ t / . To measurethe coarsening of a dry foam, we quench from an aver-age density φ = 0 . φ = − .
4. We start measuring the coarsening dynamicsafter we have stopped draining. However, our referencepoint for t = 0 is still the point of the quench. We observea scaling r ∝ t . , significantly slower than the predic-tions of von Neumann’s law (Fig 3, open circles). At latetimes in the simulation, the bubble interfaces start to be-come wet. If we restrict our analysis to times shortly aftercoarsening begins, the resulting fit exponent increases to r ∝ t . . This indicates that we are probably seeing theeffect of not having a fully dry foam at any point in time.As two-dimensional froths coarsen, they are expectedto reach a self-similar scaling state in which the normal-ized moments of the distribution of areas and of coordi-nation number are expected to become time independent.The second moment of the coordination number distri-bution has been used as a probe of this scaling state.Glazier and Weaire [11] predict that the coordinationnumber distribution should eventually reach a universallimiting scaling state with a second moment µ = 1 . µ = 2 . ± .
01 is significantly differ-ent from 1.4 (see Fig 4). In other models and simulationsof coarsening 2D foams, the observed limiting value of µ has been variously reported as 1 . . . . .
14 and 0 .
22 in magneticfroths [17].We conclude that the second moment is either not auniversal quantity, or else that there are strong transientswhich make any universal scaling regime difficult to ob-serve in practice. We note that our foam is not com-pletely dry, and that the absence of a drainage mech-anism implies that the foams becomes wetter the morethey coarsen. Whether real foams undergo a correspond-ing change of regime is not clear.We have found that in the limit of large undercooling,
FIG. 4: Evolution of the second moment of the bubble co-ordination number in a dry foam realized by the modifiedPFC equation. The bubble distribution broadens during thetransition to power-law coarsening, after which it reduces toa steady value of 2 . ± .
01. The inset shows a particulardistribution at t = 1800. the phase field crystal equation produces topologicallystabilized foams. These foams are a consequence of theresidual wavelength selection, which acts as a singularperturbation that prevents the destruction of cell wallsand forces the coarsening to take place via topological re-arrangements. The fact that this behavior emerges froma minimal model such as the phase field crystal equationsuggests a general mechanism by which foams may occurin natural systems.The ingredients of the phase field crystal equation are adriving force towards certain equilibrium densities, somesort of spatial relaxation (diffusion, viscosity, and thelike), and some competing source of wavelength selection.Such a system, in the limit where the wavelength selec-tion is weak compared to the other forces, would likelygive rise to a foam. This mathematical structure can beseen in models of magnetic froths [17], Type-I supercon-ductors [18], and in models of polygonal cells found inmelting snow [19].However, the foams produced in this limit of the PFCequation have unphysical properties. The wavelength se-lection also preserves atoms: bubbles whose diameter iscomparable to the width of the cell walls. This can beameliorated by modifying the PFC free energy to pe-nalizes atoms while encouraging stripes. As a result, thecoarsening dynamics of real foams are recovered althoughthe distribution of bubbles in the resultant scaling state appears to be somewhat different.If these differences are understood, the modified PFCequation may be a useful tool in modelling foams, as itis simpler than other existing methods such as Q-Pottsmodels [20] and minimal surface evolution [21], and isfairly easy to simulate, having only one field and a spatialstructure that is easily treated with spectral methods.We are grateful to B. Athreya for discussions aboutthe stability of the PFC equation. N. Guttenberg wassupported by a University of Illinois Distinguished Fel-lowship, and J. Dantzig acknowledges the support of theUS Dept. of Energy, under subcontract 4000076535. [1] K. Elder and M. Grant, Physical Review E , 051605(2004).[2] K. Elder, N. Provatas, J. Berry, P. Stefanovic, andM. Grant, Physical Review B , 064107 (2007).[3] N. Provatas, J. Dantzig, B. Athreya, P. Chan, P. Ste-fanovic, N. Goldenfeld, and K. Elder, JOM Journal of theMinerals, Metals and Materials Society , 83 (2007).[4] J. Mellenthin, A. Karma, and M. Plapp, Physical ReviewB , 184110 (2008).[5] P. Chan, N. Goldenfeld, and J. Dantzig, Physical ReviewE , 035701(R) (2009).[6] P. Y. Chan and N. Goldenfeld, Phys. Rev. E ,065105(R) (2009).[7] N. Goldenfeld, B. Athreya, and J. Dantzig, Physical Re-view E , 020601(R) (2005).[8] S. A. Brazovskii, Zh. Eksp. Teor. Fiz. , 175 (1975).[9] A. Bray, Advances in Physics , 481 (2002).[10] J. Von Neumann, American Society for Testing Materi-als, Cleveland (1952).[11] J. Glazier and D. Weaire, Journal of Physics: CondensedMatter , 1867 (1992).[12] F. Meyer, in (1991), vol. 2, pp. 847–857.[13] M. Marder, Physical Review A , 438 (1987).[14] H. Flyvbjerg, Physical Review E , 4037 (1993).[15] J. Glazier, M. Anderson, and G. Grest, PhilosophicalMagazine. Pt. B. Structural, Electronic, Optical andMagnetic Properties , 615 (1990).[16] J. Stavans and J. Glazier, Physical Review Letters ,1318 (1989).[17] F. Elias, C. Flament, J. Bacri, O. Cardoso, andF. Graner, Physical Review E , 3310 (1997).[18] R. Prozorov, A. Fidler, J. Hoberg, and P. Canfield, Na-ture Physics , 327 (2008).[19] T. Tiedje, K. Mitchell, B. Lau, A. Ballestad, and E. Nod-well, J. of Geophys. Res , F02015 (2006).[20] Y. Jiang and J. Glazier, Philosophical Magazine Letters , 119 (1996).[21] R. Phelan, D. Weaire, and K. Brakke, ExperimentalMathematics4