Dynamics of pattern formation in magnetic fluids

Colloids and Surfaces A. Physicochemical and Engineering Aspects, 80 (1993) 29-37 0927-7757/93/$06.00 0 1993 ~ Elsevier Science Publishers B.V. All rights reserved.

29

Dynamics of pattern formation in magnetic fluids Raymond

E. Goldstein”,*,

aDepartment bDepartment

of Physics, Joseph Henry Laboratories, Princeton University, Princeton, of Physics, Simon Fraser University, Burnaby, B.C., V5A lS6, Canada

(Received 21 August

David P. Jackson”,

1992; accepted

Stephen

A. Langerb NJ 08544,

USA

8 April 1993)

Abstract A model for the formation of labyrinthine patterns by thin domains of magnetic fluids is described. The dynamics is derived variationally from an energy functional which reflects a competition between surface tension and long-range magnetic interactions. A detailed linear stability analysis for circular shapes reveals the mechanism of pattern formation and the characteristic mode of instability. Numerical study of the fully non-linear regime reveals branched patterns whose details are sensitively dependent on initial conditions, as in experiment. Connections with glassy systems and spin models are suggested. Key words: Dynamics;

Magnetic

fluids; Pattern

formation

The existence

Introduction

of general

of these structures

questions

about

poses a number

dissipative

pattern

for-

Since the early work of Rosensweig [l] it has been known that striking patterns form when a droplet of magnetic fluid is placed in a magnetic

mation, and their remarkable resemblance to those found in three other systems (thin solid state

field. In the simplest experiment, the drop is confined between two glass plates, surrounded by an immiscible fluid such as water, and the field is applied normal to the plates. The subsequent fingering instability, an example of which is shown viewed from above in Fig. 1 (from the work of Dickstein et al. [2]), arises from a competition

films in applied

magnetic

layers of dipolar

molecules

between the surface tension at the ferrofluid-water interface and the long-range magnetic interactions. The former tends to minimize the perimeter of the domain, whereas the latter, arising from the repulsion between parallel dipoles, promotes fingering. From a conceptual point of view, this competition is further complicated by the incompressibility of the ferrofluid domain, whose area is thus constrained to remain constant during the pattern evolution. *Corresponding

author.

magnetic

domains

[3],

Type

I superconducting

fields [4,5] and monoat the air-water

inter-

face [6,7]) suggests further that there might be a common physical mechanism underlying them all. Some of the general questions which have motivated the present variational motion

study

principle

leading

configurations energetic

effects

operational?

there

quasi-two-dimensional

in

the

are kinetic

geometrical

and

these patterns?

work, we summarize

to these issues of labyrinthine

ground

To what extent

dominant

ways to describe

In the present

the interface Is the space of

and to what extent Are

Is there a

by a unique

local minima?

considerations

formation,

approach dynamics

governs

characterized

pattern

topological

which

to these patterns?

state, or degenerate are

are as follows.

a theoretical

[8,9], a model for the pattern formation in

magnetic

fluid

domains.

30

R.E. Goldstein et al./Colloids

guing and

Surfaces A: Physicochem.

connection solutions

Eng. Aspects 80 (1993) 29-37

between to

a

the observed

well-known

patterns

optimization

problem. Experimental Dickstein instability

results et al. [2] have studied

of ferrofluid

domains

the branching

in a Hele-Shaw

geometry (see also earlier work in Refs. 9 and 10). The experimental cell consists of two horizontal glass plates spaced approximately 1 mm apart by a gasket. It is filled with a surfactanttwater mixture which is immiscible with the ferrofluid drop and which wets the glass plates so that the drop can slide laterally without touching the plates. The magnetic field is produced by a large Helmholtz coil driven by a water-cooled and computer-controlled power supply. A video camera placed above the cell is used to record the dynamics of the pattern formation for later analysis. The patterns of interest are simply the ferrofluiddwater bound-

Fig. I. Experimental fingering instability of a circular ferrofluid domain after an abrupt jump of the magnetic field (top). Bottom panel shows the relaxation toward the original circular shape as the field is ramped back to zero linearly in time. Data were obtained by digitization of images collected on videotape (courtesy of Dickstein et al. 121). See also Ref. IO.

Recent experimental results [2] are first discussed, as they provide a benchmark for comparison with any theory. A simple and general variational formalism for the geometrically constrained dissipative evolution of a curve in the plane toward an energetic minimum is reviewed, and then applied to the energy functional appropriate to the magnetic fluid domains. We comment on some possible relationships between these systems and spin models in statistical mechanics, and note an intri-

aries, closed curves in the plane retrieved from the videotapes by standard edge-detection software. The shape instability of the captive drop is typically produced by ramping the field at a controllable rate to some final value. One of the easiest questions to answer with these experiments is whether the configuration space has more than one local minimum. The answer is yes, as shown in Fig. 2. Panels (a) and (b) show two distinct final states obtained from ostensibly identical initial conditions (shapes which are circular to the eye) long after the field was abruptly stepped to a low value; panels (c) and (d) show two of the many states found at a higher field. Not only are the final states of interest, but also the path which connects them to the initially circular shapes. A variational formalism The fluid motion in the Hele-Shaw geometry is viscously dominated, i.e. without inertia, and it is thus natural to study equations of motion for the ferrofluid-water interface which are first order in time derivatives. Such viscous domination also

R.E. Goldstein et al./Colloids

Surfaces

A: Physicochem.

Eng. Aspects 80 (1993) 29-37

31

viscous forces linear in the velocity, & is quadratic in I%/&. A simple model for the dissipation is one which is purely

local, i.e.

(2) where metric,

Fig. 2. Experimental stationary states of ferrofluids in the Hele-Shaw geometry long after the application of a magnetic field. The patterns in (a) and (b) are two possible final states at a given low field; those in (c) and (d) pertain to a higher field (data courtesy of Dickstein et al. [2]).

occurs in the other dipolar and magnetic systems mentioned in the Introduction, where incompressible hydrodynamics per se does not play a role. As a general dynamical model with which to investigate pattern formation with global constraints, we thus seek dynamics which dissipatively moves a shape to some minimum of an energy functional 6’. In order to address certain questions of reparameterization invariance it is useful to appeal to Lagrange’s formalism for dissipative processes [Ill]. Let CIE [0,1] be an arbitrary parameterization of the curve (e.g. arc length s), and treat the positions of the points Y(M) as the generalized coordinates qa. If 9 ({q.}, {a,}) is the Lagrangian; then the equation of motion is d a9 a2 _-_-=__ dt84, 84,

aFd 84,

(1)

The Rayleigh dissipation function & is proportional to the rate at which viscous forces dissipate energy, and its derivative with respect to the velocity 4, serves as a generalized force. For the typical

q is a friction and

the “gauge

coefficient, function”

g = Y;Y,

is the

O(cr, t) ensures

time-dependent reparameterization invariance. In the viscous limit, we neglect kinetic energy terms in the Lagrangian, and thus set 9’ = -8. Equation (1) is then rewritten in terms of functional derivatives with respect to v(a), and resealed in time to absorb 11.The result is reminiscent both of the Ginzburg-Landau model of critical dynamics [ 121 and the Rouse model of polymer dynamics [13]: (3) The incompressibility of the magnetic fluid implies that the area bounded by the curve Y(S)is constant, a constraint which may be enforced through a time-dependent Lagrange multiplier n conjugate to the area in an augmented energy functional & = &e - IIA, where G, is the energy of the unconstrained system. To determine II, it is useful to consider the general form of first-order curve dynamics dr -=uuli+wt” at

(4)

where ii and i^ are the unit normal and tangent vectors to the curve. This form serves to define U and W as the components of the velocity in the local Frenet-Serret frame. If we let the velocities U, and W, in Eq. (3) be those derived from J$ alone, then by differentiation of the augmented free energy rJ=u,+n

& we obtain w=w,+o

(5)

The constraint of area conservation can be shown to require I$dsU = 0. This determines II through a

R.E. Goldstein et al./Colloids

32

non-local

magnetic

relation

n= -;

By the reparameterization

invariance

of the shape

evolution, we may fix the as yet undetermined tangential velocity W for our convenience. A useful choice [S] is that which maintains uniform spacing on the curve.

This

relative

of the curve evolves in time as (8)

For numerical work, rather than evolving a set of vectors {Y(U)),it is often convenient to follow the tangent angles (B(a)} or the curvatures {k(a)}. These obey the equations [14] (9) and ; K(LX, t) = We now turn to the calculation functional &, and velocity function to dipolar domains.

(11)

To find &r it is useful to use an approximation of a uniform array of dipoles, thus making the problem equivalent to finding the electrostatic energy of an arbitrarily shaped parallel plate capacitor [15]. The field energy can then be written as

arc-length

gauge has W determined by U through the condition (a/&)(s/L) = 0, which yields a differential equation whose solution is

The total length

Eng. Aspects 80 (1993) 29-37

field energy 8r:

&o = yL + &f

dsU, i

of points

Surfaces A: Physicochem.

(10) of the energy U, appropriate

Energetics of dipolar domains Consider a domain in the shape of a slab of thickness h with a uniform cross-section described by an arbitrary closed curve +? in the plane. For pattern formation occurring with fixed total volume, and hence fixed area enclosed by the curve Y(S), the energy is that due to the line tension y around the boundary of the domain plus the

(12) where E is the total electric field, and the integral is over all of space. Although this form is exact, it offers little insight into the shape dependence of the energy. An important technical detail which allows the calculation to proceed is based on the work of Keller et al. [lS], in which the total field is written as E = E, + E, , where E. is the idealized uniform field which exists only between the two plates, and E, is the fringing field. Since the total field is curl free, we have V x E, = -V x E,. Noting that the curl of E, is non-zero only along the ribbon-shaped boundary of the slab, we see that E, may then be written in terms of a vector potential produced by a fictitious ribbon of current flowing around that boundary. Apart from terms proportional to the (conserved) area of the domain, &r takes the form [S] of a self-interacting boundary &t= -~~ds~du’qs).I(s.)m(R,~)

(13)

Here, Q(t) = sinh-‘(l/t) + 4 - dm. with is the distance between R = If? = Ir(s’)-r(s)1 points at arc lengths s and s’ on %?. Magnetic systems with magnetization M, are described by the substitution E, --+47cM0

(14)

In the interpretation of a self-interacting current loop, the scalar product in Eq. (13) reflects the attraction (repulsion) between parallel (antiparallel) current-carrying wires. The function @ is essentially Coulombic at large distances, but diverges only logarithmically at the origin. This cross-over

R.E. Goldstein et al./Colloids

Surfaces A: Physicochem.

Eng. Aspects 80 (1993) 29-37

33

occurs because of the finite thickness h of the slab, and prevents the integrals from diverging without additional

cut-offs.

From

these results

we deduce

that

the full energy

of Eq. (11) is determined

one

dimensionless

parameter,

number” [I] NBo = Eih*/87r*~, of the dipolar region.

i.e.

the

by

“Bond

and by the shape

0.9

Dynamics of tense interfaces Before applying the formalism to the full problem with competing surface tension and magnetic contributions, it is instructive to consider shape relaxation due to tension alone. As a physical realization of this, we can imagine an oil drop surrounded by water, prepared initially in a noncircular shape and relaxing toward a circle. It is clear that the associated energy functional &,, = yL is minimized by shapes having the smallest perimeter consistent with the prescribed area, i.e. circles. The bare normal velocity is simply U,(s) = -ylc

(15)

reflecting the usual Laplace pressure due to the curvature of the interface. The tangential force W,(s) vanishes. The Lagrange multiplier n follows from Eq. (6) as ZZ= (2n/L)y, and the dynamics reduce to two coupled differential equations for the total length and for the curvature:

0.77 0

(27.9

at

dsK-*

L

(16)

and

a+,

t)

ar

a* -Zs’+K3-$K*+gw K

(17)

with W given by Eq. (7) and z = yt. Apart from the non-locality associated with the tangential velocity, Eq. (17) is an area-conserving version of the wellknown “curve-shortening equation” [ 161, an essentially diffusive equation. It quite generally relaxes the shape to a fixed point of uniform curvature K = l/RO, with L = 271R,,, where 7~Rz is the area of

40

60

60

100

t Fig. 3. Relaxation of the length and shape (inset) of an elliptical initial condition to a circle under the area-conserving curveshortening dynamics given by Eqs. (16) and (17).

the initial relaxation

shape, as illustrated of an elliptical initial

in Fig. 3 for the condition.

Dynamics of dipolar domains Returning to the full energy functional appropriate to dipolar domains, we find the total normal velocity U,(s)=

aL -=--

20

-7x+$

P

ds’&?[dm-l]

(18) where l? = R/R is the unit vector pointing from the point s toward s’. The dipolar contribution is the Biot-Savart force due to a ribbon carrying an effective current I = Ee hc/4n around the boundary. Analytic progress can be made in the linear stability analysis about a circular shape [8,9,17,18], yielding for instance the most unstable mode in the presence of a given field. The linear analysis is of limited predictive ability regarding the final state, except with regard to the number of arms, which closely correlates with the mode number with the fastest growth. It does, however, often

R.E. Goldstein et al./Colloids Surfaces A: Physicochrm. Erg. Aspects 80 (1993) 29-37

34

give insight

into the mechanism

of the basic insta-

bility of the circular shape. The linear analysis consists

of examining

dynamics

circle, parameter-

10

of a weakly perturbed

the

ized by the polar angle cp: r(q) = C& + i(cp, Tangential

W,(cp)

velocities

(19)

also play

order, so we may choose the gauge function If the perturbation has the form i(cp, t) = 1 exp(a,t) n

zt

1

no role at linear

cos ncp

0 = 0.

(20)

<-(ferrofluids)-> 1,111

0.1

then we find that the growth

rate cn of mode n is

1

given by [S]

(monolayers)->

I II1llll

10

I 1111111 100

1000

2%/h

Fig. 4. Theoretical phase diagram dipolar domains from the linear

for pattern formation stability analysis given

of by

Eq. (21).

where $Sn(p) is a complex function of the mode number and the unperturbed aspect ratio p =

2R,lh

(22)

The aspect ratio p is very large in typical experiments, easily achieving values of 20-100 for the geometries used for a magnetic fluid, and approaching 103-lo4 for monolayers. In the limit of large p the growth rate simplifies considerably:

(- > R;

cJ,z

Y

1 -??+LV,, 2 x {(n2 - l)[l

- C - ln(n/p)]

- In n} (23)

where C z 0.577 2 15 is Euler’s constant. Grouping together terms proportional to (1 - n2) in Eq. (23), we find an effective line tension Yeff=Y

i

l-~C1-c+lnhWBo

I

(24)

which is lower than the bare tension and can be driven negative for a sufficiently large field. The logarithmic dependence on n in Eq. (23) which is destabilizing for n/p < 1, arises from the long-range interaction between the tangent vectors. A suggestive “phase diagram” in the p-N,, plane may be constructed by locating the lines along

which each of the modes n is at the maximum of CJ~.Near each of these lines we would expect the pattern to display n arms (ignoring any tip-splitting bifurcations after the initial instability). As remarked in earlier studies [ 15,191, a circular drop may become unstable either by increasing the dipolar strength or reducing the surface tension (i.e. increasing the Bond number) or by growing larger (increasing p), as shown in Fig. 4. At present, an understanding of the non-linear regime requires numerical study. Figure 5 illustrates the dynamics of fingering obtained for parameters like those in experiment. The overall structure is in good qualitative agreement with experimental patterns, displaying a well-defined finger width and the proper internal node structure. The dynamics also displays sensitive dependence on initial conditions in the sense that starting configurations in the form of perturbed circles differing only in the amplitudes of high frequency Fourier modes can lead to vastly different final states. Geometrical

and topological description of patterns

Having seen that perhaps the simplest dynamical and energetic model can reproduce the basic

R.E. Goldstein et al./Colloids Surfaces A: Physicochem. Eng. Aspects 80 (1993) 29-37

0 (4

0 (b)

Q (cl

Fig. 5. Theoretical relaxation of a circular domain to a branched final state. The material parameters are N,, = 0.54 and p = 200. Panels (a)-(f) are at equally spaced time intervals.

aspects of pattern formation seen in experiment, we now turn to some more qualitative observations regarding these shapes. Observe that each “tree” has threefold coordinated internal nodes whose angles are approximately 120”. This feature suggests the importance of surface tension in the final shapes. At the same time, the width of the fingers is quite sharply defined. In some sense, then, the patterns consist of elementary objects (fingers) linked together in any of a large number of topologically and geometrically distinct patterns. Indeed, they resemble the minimal trees which appear in the celebrated “motorway problem” [20,21] -the optimization problem of finding the connections (roads) between a given number of points (cities) such that a path exists from any given city to any other and the total road length is minimized. With the exception of a small class of spatial arrangement of the cities, the minimal path length graph is generally found by introducing additional points (called Steiner points) which are threefold coordinated nodes. Steiner minimal trees may be classified by focusing on the topology of the backbone of nodes. A

35

convenient arrangement of these trees follows [22,23] by grouping those with a given number of free ends (or equivalently, a given number of nodes) into a single generation and linking members of adjacent generations which can be obtained from one another by the addition or subtraction of a single Steiner point. Thus it is possible to organize the local minima in the configuration space of dipolar domain shapes into a hierarchical arrangement. Experimentally, it is found [2] that the degree to which this hierarchy is penetrated (i.e. the number of nodes in the final state) depends on the rate at which the magnetic field is ramped to some final value. This suggests that the nodes may be viewed as defects in the pattern created as the shape falls out of equilibrium with the applied field. The near equality found [2] in the overall scale (e.g. radius of gyration and perimeter) of patterns at a given level of the hierarchy suggests the near equality of their energies. For these reasons, the pattern formation exhibited by these magnetic fluid domains may be termed “glassy” by analogy with corresponding behavior in spin glasses [24]. A correspondence to spin systems

In the light of the glassy behavior discussed above, it is intriguing to note a correspondence between the energetics of these dipolar domains and that of a classical spin system. First, we observe that the magnetic energy in Eq. (13) is very much like an XY spin hamiltonian of the form ~ = ~ ,~. JijSi. Sj I+3

(25)

with a long-range coupling Jij between sites i and j, if we identify the tangent vector as a spin variable Si O $S)

(26)

This interpretation is not completely straightforward for two reasons. First, the arc length variables s and s’ which label two interacting tangent vectors are internal markers which do not directly specify

R.E. Goldstein et al./Colloids

36

the distance

between

the two; rather,

that distance

is the magnitude [r(s) - r(s’)J. Only for trivial shapes such as circles is there a simple relation between 1Y(S)- r(d) 1and 1s - s’l. Secondly, when comparing different final states in this light, it is not clear a priori that the number of effective spins (i.e. the perimeter) are the same, although experimentally we have found that the ensemble from

ostensibly

identical

initial

of shapes obtained conditions

has a

sharply defined perimeter [2]. Proceeding with this interpretation, we find that this spin system has vanishing net “magnetization” (not the magnetization of the domain, of course). This follows from the closure condition on the integrated tangent vector:

Surfaces

magnetic

energies.

as follows. not

valid

A: Physicochem.

Erg. Aspects 80 (1993) 29-37

Some areas for future work are

First, the local model

of dissipation

for the full hydrodynamic

problem

is of

Hele-Shaw flow, and it is thus of interest to develop a generalization of Darcy’s law for this problem (as in Ref. 18). Secondly, further study of the connection with minimal trees is clearly warranted, particularly with regard to the development of the notion

of weakly interacting

effective degrees

of freedom. Moreover, it remains to develop some analytic understanding of both the sensitive dependence on initial conditions and the kinetic control of the pattern. Finally, a more detailed understanding is needed of the connection between pattern formation in systems such as Type I superconductors and that described here.

(27) Acknowledgments the latter following from the definition of the tangent vector as 4s) = dr/ds and the condition v(L) = r(O) for a closed curve. Thus, although there is no “long-range order” in the language of spin systems, there is short-range order in the formation of arms, a process akin to the pairing of spins. It is quite remarkable that the basic feature of frustration found in many glassy systems (structural glasses, spin systems, etc.) has a clear analog here in the competition between long- and shortrange forces, leading to an apparently hierarchical arrangement of local minima in the space of configurations. This observation naturally leads to the question of whether it is possible to quantify the notion of a “metric” on the space of configurations - a distance between shapes - on the basis of the statistical mechanical definitions which have been used in the theory of spin glasses. Research is currently under way in this direction.

Discussion The results summarized here suggest that the branched patterns seen in experiment can be understood within perhaps the simplest dynamical model incorporating the competition between surface and

We thank our colleagues A. Dickstein and Erramilli for ongoing collaborations, S. and A.O. Cebers, N.P. Ong, D.M. Petrich, R. Rosensweig, M. Seul, and M.J. Shelley for important discussions. This work was funded in part by National Science Foundation Grant No. CHE91-06240, and the Natural Sciences and Engineering Research Council of Canada.

References I 2 3 4 5 6 7 8 9 10

R.E. Rosensweig, Ferrohydrodynamics, Cambridge University Press, Cambridge, 1985. A.J. Dickstein, S. Erramilli, R.E. Goldstein, D.P. Jackson and S.A. Langer, Science, 261 (1993) 1012. M. Seul, L.R. Monar, L. O’Gorman and R. Wolfe, Science, 254 (1991) 1616. F. Haenssler and L. Rinderer, Helv. Phys. Acta, 40 (1967) 659. Huebner, Flux R.P. Magnetic Structures in Superconductors, Springer-Verlag, New York, 1979. H. Miihwald, Annu. Rev. Phys. Chem., 41 (1990) 441. H.M. McConnell, Annu. Rev. Phys. Chem., 42 (1991) 171. S.A. Langer, R.E. Goldstein and D.P. Jackson, Phys. Rev. A, 46 (1992) 4894. A.O. Tsebers and M.M. Maiorov, Magnetohydrodynamics, 16 (1980) 21. A.G. Boudouvis, J.L. Puchalla and L.E. Striven, J. Colloid Interface Sci., 124 (1988) 688.

R.E. Goldstein et al./Colloids Surfaces A: Physicochem. Eng. Aspects 80 (1993) 29-37 11 12 13 14 15 16 17

H. Goldstein, Classical Mechanics, Addison-Wesley, Reading, MA, 1980. P.C. Hohenberg and B.1. Halperin, Rev. Mod. Phys., 49 (1977) 435. M. Doi and SF. Edwards, The Theory of Polymer Dynamics, Oxford University Press, New York, 1986. R.C. Brower, D.A. Kessler, J. Koplik and H. Levine, Phys. Rev. A, 29 (1984) 1335. D.J. Keller, J.P. Korb and H. McConnell, J. Phys. Chem., 91 (1987) 6417. M.E. Gage, Invent. Math., 76 (1984) 357. A.A. Thiele, Bell Syst. Tech. J., 48 (1969) 3287.

18 19 20 21 22 23 24

31

A.O. Tsebers, Magnetohydrodynamics, 23 (1988) 13; 25 (1989) 149. T.K. Vanderlick and H. Miihwald, J. Phys. Chem., 94 (1990) 886. C. Isenberg, The Science of Soap Films and Soap Bubbles, Dover, New York, 1992. E.N. Gilbert and H.O. Pollak, SIAM (Sot. Ind. Appl. Math.) J. Appl. Math., 16 (1968) 1. F.K. Hwang and J.F. Weng, Discrete Math., 62 (1986) 49. D.Z. Du, F.K. Hwang and J.F. Weng, Discrete Comput. Geom., 2 (1987) 65. K. Binder and A. Young, Rev. Mod. Phys., 58 (1986) 846.