J. Fluid Mech. (2013), vol. 714, pp. 612–633.
doi:10.1017/jfm.2012.512

c Cambridge University Press 2013

612

Magnetic field driven micro-convection in theHele-Shaw cell
K. Ērglis1 , A. Tatulcenkov1 , G. Kitenbergs1 , O. Petrichenko1 , F. G. Ergin2 ,B. B. Watz2 and A. Cēbers1 , †1 Department of Theoretical Physics, University of Latvia, Riga, Zellu 8, LV-1002, Latvia
2 Dantec Dynamics A/S, Tonsbakken 16-18, Skovlunde, 2740, Denmark(Received 2 March 2012; revised 27 September 2012; accepted 11 October 2012)

Micro-convection caused by ponderomotive forces of the self-magnetic field of amagnetic fluid in the Hele-Shaw cell under the action of a vertical homogeneous
magnetic field is studied both experimentally and numerically. It is shown that a nonpotential magnetic force at magnetic Rayleigh numbers greater than the critical valuecauses fingering at the interface between the miscible magnetic and non-magnetic
fluids. The threshold value of the magnetic Rayleigh number depends on the smearingof the interface between fluids. Fingering with its subsequent decay due to diffusionof particles significantly increases the mixing at the interface. Velocity and vorticity
fields at fingering are determined by particle image velocimetry measurements andqualitatively correspond well to the results of numerical simulations of the microconvection in the Hele-Shaw cell carried out in the Darcy approximation, which
account for ponderomotive forces of the self-magnetic field of the magnetic fluid.
Gravity plays an important role at the initial stage of the fingering observed in theexperiments.
Key words: Hele-Shaw flows, MHD and electrohydrodynamics, magnetic fluids

1. Introduction
Magnetic-field-driven micro-convection was found in the early 1980s (Maiorov &Cebers 1983) extending the previous work on magnetostatic instabilities of magnetic
liquids in the Hele-Shaw cell (Cebers & Maiorov 1980) for miscible fluids. Since theponderomotive force on the magnetizable fluid is proportional to the concentration ofmagnetic particles and the local gradient of the magnetic field strength, it is potential
only when the concentration and magnetic field strength gradients are collinear. The
field strength gradient in a homogeneous applied field results from the self-magnetic
field of the liquid. Any concentration perturbations that destroy this collinearity cause
the liquid to start to flow. An important issue that is absent in gravitational convectionis the dependence of the field gradient on the concentration field. It should be
mentioned that the equivalent of the thermal expansion coefficient, which in liquids
has an order of magnitude of 10−3 , is much larger in this case because of the strong
dependence of the magnetization of a liquid on the particle concentration.
In view of these considerations, the theoretical model of the magnetic microconvection considers the Hele-Shaw flow under the action of the ponderomotive forces† Email address for correspondence: aceb@tesla.sal.lv613

Magnetic micro-convection

due to the self-magnetic field of the fluid, the equations for the magnetostatic field,
and the diffusion equation for the concentration of the magnetic nanoparticles (Cebers
1997). The theoretical analysis shows that there is a threshold value of the criticalmagnetic field strength determined by the magnetic Rayleigh number.
Later, magnetic-field-driven micro-convection has attracted the interest of researchers
both from an experimental and a theoretical point of view. Magnetic-field-driven
fingering was observed in micro-channels using the methods of microfluidics (Derec
et al. 2008) and on the interface of a circular region concentrated by magnetic
particles (Wen, Chen & Kuan 2007). The concept of the magnetic-field-driven microconvection was experimentally studied for the case of the time evolution of gratingsof magnetic particles induced by thermophoresis at non-homogeneous illumination of
the magnetic colloids (Mezulis & Blums 2005). The theoretical analysis of the lastproblem was given by Igonin & Cebers (2003).
The theoretical analysis of the growth increments of the fingering instability inthe case of plane interface between miscible magnetic and non-magnetic fluids wasgiven by Igonin & Cebers (2003). Numerical analysis of the magnetic-field-driven
micro-convection in the frame of the described model was carried out by Chen & Wen(2002). In the work of Chen (2003) the effect of the Korteweg stress was taken into
account in the numerical simulation of magnetic micro-convection.
In spite of various efforts, the quantitative information on the main characteristicsof the micro-convection that develops is rather poor. Even the occurrence of fluid flow
during the formation of the fingers at the interface between two miscible fluids has
not yet been confirmed. Here we obtain quantitative information on the micro-flow
driven by the self-magnetic field of the magnetic colloid by applying particle image
velocimetry (PIV) methods with specially designed algorithms for processing images
with poor contrast (Ergin et al. 2010). The main obstacles for quantitative comparison
with the theoretical analysis and numerical simulations given in the first part of thework are discussed, and it is shown that in spite of small density differences between
the magnetic fluid and its solvent gravity plays an important role in the smearing of
the concentration distribution, which essentially determines the critical field strength
for the development of the magnetic micro-convection. Nevertheless, we illustrate that,in general, the numerical simulation and experimental data are in good qualitative
agreement.
2. Mathematical formulation
2.1. Model

We consider two miscible fluids where the first is a magnetic fluid and the second isa simple non-magnetic fluid. Fluids are confined in a horizontal Hele-Shaw cell anda magnetic field is applied perpendicular to the cell, the sketch of the cell is shownin figure 1. The viscosities of the fluids are equal. The ponderomotive forces of the
non-homogeneous self-magnetic field on the interface between fluids cause fingering
instability. Its growth is described by the set of equations, which includes the Darcy
equation, the convection–diffusion equation and equations for the magnetostatic field(Cebers 1997; Igonin & Cebers 2003) and reads
−∇p −

12η
2M(c)
u−
∇ψm (c) = 0,h2h∂c
+ (u · ∇)c = D∇ 2 c,∂t

∇ · u = 0,

(2.1)
(2.2)

614

K. Ērglis and others
H
Ferro liquidhLiquidzy0x

F IGURE 1. The Hele-Shaw cell sketch.

where p is pressure, u = (ux (x, y), uy (x, y)) is the depth averaged velocity, η is theviscosity of the fluid, h is the thickness of the Hele-Shaw cell, D is isotropic constantdiffusion coefficient and c is the concentration of magnetic fluid normalized by its
value far from the interface. The magnetization M(c) is taken to be proportional to theconcentration of the magnetic fluid c (M = M0 c) and the value of the magnetostatic
potential ψm on the boundary of the Hele-Shaw cell is given by Cebers (1981) andJackson, Goldstein & Cebers (1994)
Zψm (r, t) = M0 c(r0 , t)K(r − r0 , h) dS0 ,(2.3)
where the integration
p is performed over the boundary of the Hele-Shaw cell,K(r, h) = 1/|r| − 1/ |r|2 +h2 .
Equations (2.1) and (2.2) extend in the Darcy approximation the model which isusually applied for the study of the displacements of miscible, non-magnetic fluids inthe Hele-Shaw cell (Goyal & Meiburg 2006).
The boundary conditions for the velocity components and the concentration of thefluid are as follows:ux (0, y) = uy (0, y) = 0, c(0, y) = 1,
ux (Lx , y) = uy (Lx , y) = 0, c(Lx , y) = 0,

(2.4a)(2.4b)and the conditions of the periodicity across the Hele-Shaw cell areu(x, 0, t) = u(x, Ly , t),c(x, 0, t) = c(x, Ly , t).

(2.5)

The boundary conditions (2.4) require that the fluid is motionless at both endsof the cell. The motion of the liquid arises from a non-potential magnetic force−2M(c)∇ψm /h.
The equations are put in dimensionless form by introducing the following scales:
length h, time h2 /D, velocity D/h and magnetostatic potential M0 h. As a result the setof dimensionless equations reads−∇p − u − 2Ram c∇ψm = 0, ∇ · u = 0,∂c
+ (u · ∇)c = ∇ 2 c.
∂t

(2.6)
(2.7)

615

Magnetic micro-convection(a) 0.30(b) 0.30

0.25

0.25

0.20

0.20

0.15

0.15
0.10

0.10

0.05

0.05

0
0

5

10

15

20

0

5

10

15

20

k

k

F IGURE 2. (a) The neutral curves of the magnetic micro-convection for different values ofthe smearing parameter (t0 = 0.5, t0 = 0.1, t0 = 0.01, t0 = 0.0001). (b) The theoretical neutral
curve (solid line) and the neutral curve calculated numerically (dashed line).

Here Ram = M02 h2 /12ηD is the magnetic Rayleigh number determined by the ratio of
the characteristic time of the diffusion τD = h2 /D and the characteristic time of motion
due to non-homogeneous self-magnetic field of the fluid τM = 12η/M02 .
2.1.1. The linear stability analysis
In order to test the numerical algorithm, we compare the theoretical results of alinear stability analysis with the numerical results. The solution of the quasi-steadylinear stability problem for miscible magnetic fluid in the Hele-Shaw cell is foundin Igonin & Cebers (2003). An analytical solution may be found for the value ofthe smearing parameter t0 = 0. The linear perturbation of a quiescent base stateis represented by {vx , vy , c, ψm }(x, y, t) = {0, 0, c0 , ψm0 }(x) + {vx0 , vy0 , c0 , ψm0 }(x)eiky+λt√(c0 = 0.5(1 − erf(x/2 t0 ))), where the parameter t0 characterizes the thickness of
the interfacial region. We focus on the linear analysis of the stability of the interface inthe formal limit t0 = 0 when the concentration distribution is step-like. The dispersion
relation reads

ks + Ram [2J(s, k) − kf (k(s + 1))] = 0,(2.8)
pwhere the parameter s is s = 1 + λ/k2 and the functions J(p, q) and f (p) are definedby the integrals
Z ∞pJ(p, q) =e−pz (K0 (z) − K0 ( z2 + q2 )) dz,(2.9)
0f (p) =

Z ∞0e−pz ln(1 + z−2 ) dz.

(2.10)

Here K0 is the modified Bessel function of the second kind (McDonald function).
Neutral curves shown in figure 2(a) are calculated numerically as follows.
The Fourier coefficients ĉn for the concentration field averaged along the x-axisRPNy
hc0 i(y) = Lx c0 (x, y) dx/Lx = n=0hĉn (t)ieikn y are fitted by hĉn i = hĉ0n ieλn t , and thegrowth increments λn are found in dependence on the wavenumber kn . The neutral

616

K. Ērglis and others
0.24
0.22
0.20
0.18
0.16
0.14
0.12
0.10
0.08
0.06

10 –4

10 –3

10 –2

10 –1

10 0

F IGURE 3. The inverse values of the critical magnetic Rayleigh number in dependence on thesmearing parameter t0 . Both triangles shown are analytical results from Igonin (2004), and theblack circle shows analytical results for the sharp front.

curves for definite values of the smearing parameter t0 are found by solving theequation λn (Racm , kn ) = 0 by a spline interpolation.
To reproduce in the simulations growth increments theoretically calculated fora sharp diffusion front it is necessary to take the concentration perturbation
corresponding to the theoretically found eigenfunction. The critical values of themagnetic Rayleigh number in dependence on the wavenumber for t0 = 0.0001 with
perturbation corresponding to the eigenfunction for a sharp front are shown infigure 2(b) by the dashed line. They agree reasonably well with the theoretical values
obtained from (2.8) for t0 = 0 and shown in figure 2(b) by the solid line. Their small
difference is due to difference of eigenfunctions for sharp and smeared fronts.
The inverse values of the critical magnetic Rayleigh number obtained numericallyare shown in figure 3 in semi-logarithmic coordinates as a function of t0 . They show
agreement with the available theoretical results (Igonin 2004) for sharp and smeared
fronts. These tests validate the numerical algorithm.
The numerical algorithm described in appendix A is applied further to study
fingering at magnetic micro-convection in the nonlinear regime.
2.2. Numerical results
The initial concentration perturbation at the interface is introduced by applying relation(A 6) with ς(y) taken to be random. In long-time numerical calculations typically
512 × 512 collocation points are used for the space discretization. The time step 1ttypically is 1t = 10−5 .
Typical fingering patterns for three values of Ram and the smearing parametert0 = 0.05 are shown in figure 4. The pattern formation has several stages as may beseen in figure 4. At first (t = 0.0014 and t = 0.002) well-resolved fingers are formedthat later coarsen with time (t = 0.0035 and t = 0.009). Finally, the finger patternsmears out because of diffusion of the magnetic particles (t = 0.013).

Magnetic micro-convection

617

F IGURE 4. Concentration images for magnetic Rayleigh number Ram = 300, Ram = 500,
Ram = 1500, at time t = 0.0014, t = 0.002, t = 0.0035, t = 0.009, t = 0.013, with Lx = 6,
Ly = 6 and t0 = 0.05.

Snapshots of the concentration, vorticity, stream function and velocity fields for
t0 = 0.05 and Ram = 500 at time moment t = 0.013 are shown in figure 5. The
velocity field is visualized by using the line integral convolution technique. Line
integral convolution (LIC) is a powerful technique for generating striking images
from vector data and for visualizing vector fields. Introduced by Cabral & Leedom(1993), the method has rapidly found many application areas, ranging from computer
art to scientific visualization. Given a vector field and an input texture, line integral
convolution produces an output texture in which the data values are highly correlated
in the direction of the flow. LIC produces an intensity output image, which is
reproduced as a greyscale picture. This image includes tangential information of vector
field as well as information on its magnitude. The formation of alternating array ofvortices may be seen at figure 5(b,c).
To make a qualitative comparison with the experimental data, described in the nextpart, the time dependence of the maximal vorticity for several values of the magnetic
Rayleigh number is calculated and shown in figure 6. The characteristic maximum
corresponds to the transition from the development of the magnetic micro-convection
to its decay due to particle diffusion. The spiky character of the curve at larger valuesof the magnetic Rayleigh number is due to the emerging vortices, which intensify themixing.
3. Experimental demonstration
Magnetic micro-convection is studied on the interface between two miscible
fluids. The first fluid is distilled water which contains 0.1 % fluorescent tracer
micro-particles with diameter 1 µm for PIV measurements. The second is a waterbased magnetic fluid with γ –Fe2 O3 particles of average radius 14 nm, saturation

618

K. Ērglis and others

(a) 6

0.9

(b) 6

5

4000

5
2000

0.7
4

4
0.5

3

0

3

2

2

–2000

0.3
1

1
0.1

0

1

2

3

4

5

0

6

–4000
0

1

2

3

4

5

(d)

(c) 6

60

15

5
4

10

50

5

40

3

0

2

–5

30
20

–10
1

10

–15
0

6

1

2

3

4

5

6

F IGURE 5. Snapshots of concentration (a), vorticity (b), stream function (c) and velocity
field (d) for miscible magnetic and non-magnetic fluids in the Hele-Shaw cell. Snapshotscorrespond to the magnetic Rayleigh number Ram = 500 at time t = 0.013, with Lx = 6,
Ly = 6 and t0 = 0.05.
11
10
9
8
7
6
5
4
3
–7

–6

–5

–4

–3

–2

–1

log(t)F IGURE 6. The relaxation of the maximal value of vorticity field ωmax for different values ofmagnetic numbers Ram on a logarithmic scale.

magnetization 10 G and initial susceptibility 0.064, as determined by the vibrating
sample magnetometer measurements. The density of the liquid is % = 1.077 g cm−3 .
The diffusion coefficient of the particles is determined by diffusion light scattering andequals D = 5.5 × 10−8 cm2 s−1 . Magnetic nanoparticles of γ –Fe2 O3 are synthesizedMagnetic micro-convection

619

F IGURE 7. Determination of the threshold value of the magnetic field for finger-like patterndevelopment.

by a precipitation reaction (Massart). The Hele-Shaw cell is prepared from two
microscope cover glasses separated by a 127 µm thick Parafilm. The Parafilm is
cut to form channels for air and fluids. The horizontal dimensions of the fluid channel
are 18 mm × 6 mm. After assembly the cell is heated to weld the cover glasses. A
few drops of magnetic fluid are placed on one side of the cell and one half of thecell is filled with magnetic fluid because of capillary forces. Then the magnetic field is
switched on, and a small amount of the non-magnetic fluid is added to the other sideof the cell using a small tube and a syringe with a microscrew until the solvent comes
into contact with the magnetic fluid.
The homogeneous magnetic field with strength H ∈ (0–300) Oe is created in thecentral part of the coil where the Hele-Shaw cell is placed. The coil with diameter
34 mm and height 16 mm consists of 38 turns of 0.7 mm wire.
The theoretical analysis (Cebers 1997; Igonin & Cebers 2003) and numericalcalculations described in Part 2 show that the smearing of the diffusion front between
the two miscible fluids plays an important role in the development of the magnetic
micro-convection.
The experimental data shown in figure 7 clearly confirm the existence of the criticalvalue of the field strength below which the fingering due to the magnetic microconvection is not observed. If we estimate the critical field strength according to thedata shown in figure 7 to be 3.4 G, the magnetic Rayleigh number for the fluid withthe given physical properties is Ram = 1150, a value much larger than the critical value
corresponding to t0 = 0 and shown in figure 3. Figure 3 shows that such large valuesof the critical magnetic Rayleigh number correspond to large values of the smearing

620

K. Ērglis and others

0

0.2

0.4

0.6

0.8

1.0

F IGURE 8. Concentration field at t = 0.1, 1, 5 and Ra = 400.

parameter, which should not be expected if the particle diffusion were dominant. Largesmearing may be explained by the effect of the gravitational force. Although thedensities of the magnetic and carrier liquids are similar, nevertheless a simple estimateshows that the gravitational Rayleigh number Ra = 1%gh3 /8Dη at a density differenceof 0.077 g cm−3 is quite large and exceeds 3.5 × 104 . As a result, after the liquids
make contact counterflows arise in the lower and upper parts of the Hele-Shaw cell.
Their role in the smearing of the interface may be illustrated by the Stokes modelwith the concentration c dependent gravity force and the diffusion equation. Scaling
time by h2 /4D, length by h/2 and the velocity by 1ρgh2 /4η (h is the thickness of theHele-Shaw cell), the dimensionless set of partial differential equations (PDEs) reads−∇p + ∇ 2 v − cey = 0,
(3.1)∂c
+ Ra(v · ∇)c = ∇ 2 c.
(3.2)∂t
The evolution of the concentration field and the liquid flow are calculated by using
the Comsol software with the initial concentration distribution given by c0 (x, t0 ) at
t0 = 0.025:!
Z x/(2√t0 )1
2
c0 (x, 0) =1− √exp(−y2 ) dy .
(3.3)
2
π 0
The results for Ra = 400 (the results for other values of the Rayleigh number are
qualitatively the same) in the cell x ∈ [−10 : 10] and y ∈ [−1 : 1] are shown infigures 8 and 9. The calculated concentration distribution in the middle of the cell(y = 0) at different time moments and the evolution of the horizontal component of
the velocity in the middle cross-section of the cell (x = 0) are shown in figures 10
and 11, respectively. We see that the flow due to gravity smears out the concentration
distribution and, after some transition time, the concentration distribution with a muchwider transition layer is established, as may be seen in figure 8. The time dependenceof the maximal value of the velocity of the arising flow for different values of the621

Magnetic micro-convection–0.020–0.015–0.010–0.005

0

0.005

0.010

0.015

0.020

F IGURE 9. Evolution of the horizontal velocity field at t = 0.1, 1, 5 and Ra = 400.
1.0

0.8

0.6

c
0.4

0.2

0
–10

–8

–6

–4

–2

0

2

4

6

8

10

x

F IGURE 10. Smearing of the concentration distribution due to the flow caused by gravity at
t = 0.1, 1, 5 and Ra = 400.

Rayleigh number is shown in figure 12. Thus, we see that gravity causes considerable
smearing of the interface before the fingering due to the magnetic micro-convection
starts. This result explains the rather large values of the magnetic Rayleigh number
necessary for its initiation.
Estimating the characteristic velocity of the flow caused by gravity according
the data in figure 12 by v = 10−2 1%gh2 /4η in dimensional units that give thecharacteristic velocity of flow as 31 µm s−1 , we see that it is close to the values
determined from the contrast variation of the interface observed in the experiment andshown in figure 13. The characteristic times of the decay of the flow caused by gravity
0.1h2 /4D ' 70 s are also in agreement.

622

K. Ērglis and others
0.025

0.015

0.005

–0.005–0.015–0.025
–1.0

–0.5

0

0.5

1.0

y

F IGURE 11. Horizontal component of velocity in the middle cross-section of the cell fordifferent dimensionless time moments t = 0.1, 1, 5 and Ra = 400.
0.04

0.03

0.02

0.01

0

2

4

6

8

10

t

F IGURE 12. Time dependence of the maximal horizontal velocity withRa = 400, 1000, 5000.

The experimental results of the development of magnetic micro-convection forseveral values of the magnetic field strength are shown in figure 14. The patterns
are similar to those observed in numerical experiments for moderate values of themagnetic Rayleigh number and shown in figure 4. The grainy structure of the interface,which may be seen in figure 14 for larger values of the field (18 and 27 G), iscaused by the normal field instability (Rosensweig 1985) of the layer of magnetic fluid
formed by the gravitational convection.
The estimates that confirm that the grainy structure, which may be seen in figure 14
for larger values of the field is caused by the normal field instability are given inappendix B. The results for the critical field strength and the structure period given
the physical parameters of liquids described above are in good agreement with the623

Magnetic micro-convection
700.08600.08x (µm)

500.08400.08300.08200.08100.08000.08

0

10

20

30

40

50

60

70

80

t (s)

F IGURE 13. (Colour online) Experimental data of the propagation of the interface at amagnetic field value of 2.25 G.

F IGURE 14. Fingering instability development in dependence on the magnetic field strength.

experiment (figure 14). In order to check that the stratification resulting from the flow
due to gravity and subsequent normal field instability indeed cause the grainy structure,a glycerol–water mixture with density matched to the density of the magnetic fluid
was prepared. The fingering at the interface of magnetic and non-magnetic fluids withmatched densities for field of 27 G is shown in figure 15. We clearly see that thegrainy structure is absent now. This result and the estimates given in appendix B
confirm that the grainy structure is caused by the normal field instability on theinterface between two miscible fluids. It should be mentioned here that matching
the density of the magnetic fluid is a subtle issue since it is impossible to use
salt solution which causes the coagulation of the electrostatically stabilized magnetic
fluid. Using the glycerol–water solution creates another problem: the viscosity of theglycerol–water solution is higher than the viscosity of the ferrofluid. Nevertheless,

624

K. Ērglis and others

F IGURE 15. Fingering at the interface between the fluids with matched densities, with
B = 27 G.

F IGURE 16. Vector field at time 8.5 s and B = 18 G.

we believe that the reason of the grainy structure formation now is confirmed boththeoretically and experimentally. Some small-scale stripe-like structure which may beseen in figure 15 is caused by the fingering phenomenon on the borders of theHele-Shaw cell wetted by the ferrofluid.
Thus, the full model of the magnetic micro-convection is much more complicated asconsidered here and should incorporate three-dimensional equations for the magnetic
field, concentration and flow.
Quantitative information about the velocity field of the micro-convection is obtained
by PIV measurements. The main problem of using PIV algorithms for magnetic microconvection is that it is impossible to determine the fluid flow in the darker regions
because of poor image contrast. Therefore, special image-processing algorithms
had been developed to improve the contrast (Ergin et al. 2010). This algorithm
includes image preprocessing for local image normalization and a Gaussian filter
in order to reduce noise. The vector field after preprocessing is obtained by anadaptive correlation algorithm with 2 refined steps, 2 passes per step, 32-pixel-wide
interrogation windows with 50 % overlap. The velocity field of the fingering patternobtained at H = 180 Oe and t = 8.5 s is shown in figure 16. Vortices formation similar
to that observed in numerical experiments (figure 5) is visible.
The values of the fingering velocity for one value of the field strength are
shown in figure 17. Maximal value of the fingering velocity may be compared625

Magnetic micro-convection
140
120
100
80
60
40
20
0

2

4

6

8

10

12

14

16

18

20

t (s)

F IGURE 17. (Colour online) Finger velocity in dependence on time at a magnetic field valueof 6.75 G.
1200
1200
1000
800
800
600

400
0

500

1000

1500

400
200
0
10–3

10–2

F IGURE 18. Finger velocity in dependence on time at different values of the magnetic
Rayleigh number Ram = 300, 500, 1000 and 1500, with t0 = 0.05.

with data obtained from numerical calculations. Finger velocities in dependence ontime obtained from the numerical simulations at several values of the magnetic
Rayleigh number are shown in figure 18. Its maximal value may be fitted bymax(vx ) = 0.72Ram . As a result, the dimensional velocity of fingering for thephysical parameters of the ferrofluid given above at field strength 6.75 G is equal
to 142 µm s−1 : a value very close to that obtained in the experiments (figure 17).
By using the obtained velocity field, vorticity values are calculated with the Dantec
dynamics software, which uses the neighbouring window velocity value difference
to find the vorticity for each window. Afterwards we export the vorticity fieldvalues and post process them using MATLAB to extract the necessary information,e.g. the maximum absolute vorticity time dependence, shown in figure 19. Thus, thedeveloped algorithm allows us to obtain the evolution of the vorticity field during thedevelopment of the magnetic micro-convection and its decay due to the diffusion of
the magnetic nanoparticles. It may be characterized by a rapid increase of the vorticity

626

K. Ērglis and others
G
G
G
G
G
G

101

100

10 –1

100

101

102

F IGURE 19. Vorticity dependence on time for different magnetic field values.

–0.3

–0.4

–0.5

–0.6

–0.7

0

5

10

15

20

25

30

H (G)

F IGURE 20. The exponent of the power law α for different magnetic field values.

followed by its slow decay. Again a qualitative similarity between the numerical
data figure 6 and experimental data figure 19 is observed. Experimental data for thevorticity decay may be fitted by the power law ωmax ∼ tα . The values of the exponent
α obtained are shown in figure 20. This figure illustrates the increase of the exponent
with the field strength, which corresponds to the slower decay of the fingering atlarger fields. This trend correlates with the increase of the exponent obtained in
the numerical simulations at several values of the magnetic field strength shown infigure 6.
4. Conclusions
The proposed model of the magnetic micro-convection qualitatively describes theexperimental data on the development and decay of the magnetic micro-convection
as shown by the results of numerical simulations. The characteristic feature of themagnetic-field-driven micro-convection is an initial fast development of the fingeringMagnetic micro-convection

627

with a subsequent decay of the flow due to the diffusion of the nanoparticles. As aresult considerable enhancement of the mixing is achieved, which may be interestingfor some applications of the magnetic nanoparticles in microfluidics. The smearing of
the diffusion front caused by the action of gravity force because of the small density
difference between the magnetic colloid and its solvent plays an important role in thethreshold value of the magnetic field strength for the development of the magnetic
micro-convection. Stratification caused by gravitational convection is the reason for thedevelopment of the normal field instability. The theoretical analysis of the normal fieldinstability on the interface between two miscible fluids shows good agreement with theexperimental observations.
Acknowledgement
This work was supported by the ESF project 2009/0223/1DP/1.1.1.2.0/APIA/
VIAA/008.
Appendix A. Numerical algorithm
In order to apply spectral methods and the fast Fourier transform (FFT) we represent
the concentration by a sum of two parts:

c(x, y, t) = cD (x, t) + c0 (x, y, t),(A 1)

ψm (x, y, t) = ψmD (x, t) + ψm0 (x, y, t),(A 2)

where cD (x, t) is given by the solution√ of the one-dimensional convection–diffusion
equation as cD (x, t) = 0.5(1 − erf(x/ 4t)) and corresponds to a step-like concentration
at t = 0. The concentration perturbation c0 is a periodic function in the x, y directions(c0 (x, y) = c0 (x + Lx , y + Ly )). As a result, the numerical solution should be found forthe periodic function c0 (x, y).
Similarly, the magnetostatic potential is split into two contributions
hereψmD (x, t) =
ψm0 (r) =

Z +∞
−∞

Z

S

cD (x − ζ, t) ln(1 + ζ −2 ) dζ,c0 (r0 , t)K(r − r0 , 1) dS0 .

(A 3a)(A 3b)

Integration in (A 3a) is carried out through all of the infinite region occupied by themagnetic fluid. The fast decay of the field strength ∂ψmD /∂x with the distance from thediffusion front justify the use of periodic boundary conditions. Integration in (A 3b) is
carried out along the cell with dimensions Lx and Ly .
The Fourier component of ∂ψmD /∂x is obtained by using the convolution theorem
according to the relation:∂ ψ̂mD∂ ĉD
1 − e−|kn |(kn ) = 2π
(kn )
,
∂x
∂x
|kn |

(A 4)

where the Fourier coefficients ∂ ĉD /∂x are calculated by a one-dimensional FFT (thedomain size is increased in the x direction from Lx to 4Lx to eliminate the Gibbs
phenomenon).
The perturbation of the concentration at the initial time is introduced by asmall displacement ς (y) of the isoline c = 0.5 at x = x0 (for all numerical

628
K. Ērglis and others
calculations x0 = Lx /2) according to the relation
x − x0 + ς(y)1√1 − erfc(x, y, t0 ) =
,
2
4t0

(A 5)

where t0 is the parameter that characterizes the smearing.
P
The small perturbation ς is defined by the Fourier series ς(y) = Nn=1 [an cos(qyn) +bn sin(qyn)] (q = 2π/Ly ).
For the test of the results of linear analysis
1x − x0c(x, y, t0 ) =1 − erf √+ ς(x, y),
(A 6)
2
4t0
where the amplitudes of perturbation modes an , √bn are defined by eigenfunctions−qn 1+λ/(n2 q2 )|x−x0 |(λ is given by thefound theoretically and read {an , bn } = {εa , εb }edispersion relation (2.8) for fixed value of magnetic Rayleigh number Ram ).
Equations (2.6) and (2.7) are solved numerically by the spectral method in thevorticity stream function formulation (Tan & Homsy 1988; Zimmerman & Homsy
1992). The stream function ψ is defined as ux = ∂ψ/∂y and uy = −∂ψ/∂x, and thevorticity ω is
ω = −∇ 2 ψ.

(A 7)

Taking into account (A 1) the convection–diffusion equation (2.7) can be written as
follows∂c0 ∂ 2 c0 ∂ 2 c0 ∂ψ ∂cD ∂c0∂ψ ∂c0+= 2 + 2 −
+
.
(A 8)
∂t
∂x
∂y
∂y
∂x
∂x
∂x ∂y
Equation (2.6) gives the vorticity equation∂ψm ∂c∂ψm ∂cD ∂c0+−.
ω = −2Ram
∂y
∂x
∂x
∂x ∂y

(A 9)

In order to solve numerically (A 8) and (A 9) the Fourier spectral method is used.
The concentration perturbation c0 , stream function ψ and vorticity ω are presented by
the Fourier seriesc0 (x, y, t) =ψ(x, y, t) =ω(x, y, t) =

N−1 X
M−1
X

ĉnm (t)ei(kn x+qm y) ,(A 10a)ψ̂nm (t)ei(kn x+qm y) ,(A 10b)ω̂nm (t)ei(kn x+qm y)(A 10c)n=0 m=0
N−1 X
M−1X
n=0 m=0
N−1 X
M−1X
n=0 m=0

where kn = 2πn/Lx , qm = 2πm/Ly are wavenumbers (n, m = 0, 1, 2 . . .). Functionsĉnm , ψ̂nm and ω̂nm are calculated by using FFT in collocation points xn = (Lx /N)nand ym = (Ly /N)m (n = 0, 1, 2 . . . N − 1, m = 0, 1, 2 . . . M − 1). The total numberof collocation points N, M is proportional to 2p , p = 1, 2, 3 . . . and depend on theconcrete case.

629

Magnetic micro-convection

The nonlinear terms in (A 8) and (A 9) are given by applying a FFT:
N−1 M−1∂ψ ∂cD ∂c0∂ψ ∂c0 X X
J(x, y, t) =Ĵnm (t)ei(kn x+qm y) , (A 11a)
+
−
=
∂y
∂x
∂x
∂x ∂y
n=0 m=0


N−1 X
M−1
0
0
X
∂ψm ∂cD ∂c∂ψm ∂cQ(x, y, t) =+−=Q̂nm (t)ei(kn x+qm y) . (A 11b)
∂y
∂x
∂x
∂x ∂y
n=0 m=0
As a result, the equations for the amplitudes of the Fourier modes readω̂nm = (kn2 + q2m )ψ̂nm ,(A 12)

ω̂nm = −2Ram Q̂nm ,(A 13)∂ ĉnm= −(kn2 + q2m )ĉnm − Ĵnm .
(A 14)∂t
The values of the stream function ψ(tj ) for given values of the concentration c(tj )are found from (A 12) and (A 13). Since ψm0 is periodic in both the x and y directions,ψ̂m0 (k) is obtained as follows:ψ̂m0 (k) = 2πĉnm (k)T̂nm (k),(A 15)

where k = (kn , qm ) is the two-dimensional wave vector and |k| = kn2 + q2m . The
Fourier coefficients ĉnm are calculated using a two-dimensional FFT algorithm andT̂nm (k) = (1 − e−|k| )/|k| is the Fourier transform of K(r, 1). All derivatives of ψm0 are
calculated spectrally from the obtained values of ψ̂m0 (k).
The linear equation (A 14) is solved for the known stream function by applying
the linear propagator method, which factors out the leading-order linear term prior to
2 2
2 2
the discretization. Introducing χ̂nm = e(kn +qm )t ĉnm the equation ∂ χ̂nm /∂t = −e(kn +km )t Ĵnm
is discretized by using the three-step Adams–Bashforth method (Samarskij & Gulin1989). The result readspĉ(j+1)∗= ĉjnm γ −nm1t2(j−2) 3(23Ĵjnm γ − 16Ĵ(j−1)
nm γ + 5Ĵnm γ ),122 2
γ = e−(kn +km )1t .

(A 16)(A 17)(j+1)∗(j+1)∗
Further, corrected values of ψnm
and Ĵnm
are found from (A 12), (A 13) and(A 11a), respectively. As a result, the Fourier components of the concentration are
obtained for the next time step

1t 2
1t (j+1)∗(Ĵnm + Ĵjnm ) −(k + km2 )(ĉ(j+1)∗+ ĉjnm ).
(A 18)nm
22 n
The concentration field is obtained by applying the inverse FFT to ĉ(j+1)
and addingnmcD (t). Since the numerical scheme used is not fully implicit in time, the solution may
show numerical instability.
jĉ(j+1)
nm = ĉnm −Appendix B. Normal field instability at the interface of miscible magnetic andnon-magnetic fluids
Normal field instability of the free horizontal surface of a magnetic liquid is wellunderstood. However, this instability has not been considered yet in the case of theinterface between two miscible fluids. In this appendix it is shown that in the case

630
K. Ērglis and others
when t0 → 0 it is possible to obtain analytically the dispersion equation for the growthincrement of the instability.
The Stokes equation in the case of miscible liquids and the diffusion equation forthe particle concentration read
 2H
− 1ρcgez = 0, ∇ · v = 0,
(B 1)−∇p + η∇ 2 v + χ0 c∇
2∂c+ v · ∇c = D∇ 2 c,
(B 2)∂twhere the z-axis is along the vertical direction, c is the concentration of magnetic
particles normalized to 1 by the concentration far from the interface, and 1ρ is thedensity difference of the magnetic and non-magnetic liquids.
At equilibrium when the magnetic field is parallel to the gradient of theconcentration c0 it reads
H = H0 = (0, 0, H∞ /(1 + 4πχ0 c0 )).

(B 3)Here H∞ is the magnetic field strength far from the interface in the non-magnetic
liquid and χ0 is the magnetic susceptibility.
The equation for the magnetic field in the case of small perturbations aroundc = c0 (z) reads (with c = c0 + c0 , H = H0 + H 0 and H 0 = ∇ψ 0 )∂∂∂((1 + 4πχ0 c0 )Hx0 ) + (4πχ0 c0 H0z ) + ((1 + 4πχ0 c0 )Hz0 ) = 0.
∂x∂z∂z
In the limit when c0 is step-like, it follows from (B 4) thatH = (1 + 4πχ0 c0 )Hz0 +4πχ0 c0 H∞
1 + 4πχ0 c0

(B 4)

(B 5)

is continuous on the interface where c0 (0+ ) = 0 and c0 (0− ) = 1.
The equations for the velocity and concentration perturbations readand− ∇p0 + η∇ 2 v + χ0 c0 ∇(H0 · H 0 ) − 1ρc0 gez = 0

∂c0
dc0
+ vz
= D∇ 2 c0 .
∂tdz
Looking for perturbations periodic in the horizontal direction
we have(v0 , c0 , ψ 0 ) = (v(z), c(z), ψ(z)) exp(ikx),
 2

∂c
dc0
dc
2
+ vz
=D
−
k
c
,
∂t
dz
dz2
 2
2
χ0 k2 H∞ dc0
d
2
η
−
k
v
=
−
H − 1ρcgk2 ,zdz2
1 + 4πχ0 c0 dz


d
dψ
d
(1 + 4πχ0 c0 )− (1 + 4πχ0 c0 )k2 ψ + (4πχ0 cH0z ) = 0.
dz
dz
dz

(B 6)

(B 7)

(B 8)

(B 9)
(B 10)(B 11)Rδ
For step-like c0 (z) the boundary conditions are obtained by integration −δ dz(· · ·)of (B 9)–(B 11) across the diffusion layer and taking the limit δ → 0. They

Magnetic micro-convection

631

read (µ∞ = 1 + 4πχ0 )c(0+ ) = c(0− ),




dc +
dc −



D (0 ) − D (0 ) = −vz (0),


dz
dz



+
−

vz (0 ) = vz (0 ),





dvz +
dvz −


(0 ) =
(0 ),



dz
dz
2
2
d vz +d vz
(0 ) = 2 (0− ),


dz2
dz



3
3
2

d vz +
d vz −(µ∞ − 1)k H∞


(0
)
−
(0
)
=
H(0),



dz3
dz3
4πηµ∞



+
−

ψ(0 ) = ψ(0 ),




dψz +
dψz −
(µ∞ − 1)c(0)H∞ (0 ) + (µ∞ − 1)c(0)H∞ = µ∞(0 ) +.dz
dzµ∞(B 12)

Equations (B 9)–(B 11) at boundary conditions (B 12) are analysed in the quasistationary case when v, c, ψ ∼ exp(λt). The solutions of (B 9)–(B 11) that satisfy thecondition v, c, ψ → 0 at |z| → ∞ and the boundary conditions for the magnetostatic
potential and concentration readpc = E exp(− λ/D + k2 z) (z > 0),



p


2

c = E exp( λ/D + k z) (z < 0),





(µ∞ − 1)H∞ E p
22( λ/D + k (µ∞ + 1) + (µ∞ − 1)k) exp(−kz)ψ =−λ/Dµ∞ (µ∞ + 1)



p

(µ∞ − 1)H∞ p


2
2
+
λ/D + k exp(− λ/D + k z)E (z > 0),



λ/D


(µ∞ − 1)H∞ E p
2
2(B 13)
ψ=( λ/D + k (µ∞ + 1) − (µ∞ − 1)k) exp(kz)λ/Dµ∞ (µ∞ + 1)




p
(µ∞ − 1)H∞ p


2
2

−
λ/D
+
k
exp(
λ/D
+
k
z)E
(z
>
0),

2

λ/Dµ∞



2

p

1ρgk E

2 z)

exp(−
λ/D
+
k
(z
>
0),
vz = C exp(−kz) + Dz exp(−kz) −2η (λ/D)




2
p

1ρgk E2
vz = A exp(kz) + Bz exp(kz) −
exp(
λ/D
+
k
z)
(z
<
0).
2
η (λ/D)
The remaining five boundary conditions give the set of five equations for the
unknown set of constants A, B, C, D, E. Its solubility condition gives the dispersion
equation for the growth increment of small perturbationspp
2
4 λ/D + k2 (µ∞ − 1)2 (µ2∞ + 1)H∞λ/D + k2 − k−k
8πηDµ2∞ (µ∞ + 1)
λ/Dk
p
2
3
1ρg λ/D + k2 (2k − λ/D) − 2k=.
(B 14)Dηk (λ/D)2

632

K. Ērglis and others
3.5
3.0
2.5
2.0
1.5
1.0
0.5
0

0.2

0.4

0.6

0.8

1.0

1.2

1.4

1.6

k

F IGURE 21. Growth increments of the normal field instability at the interface of two miscible
fluids in dependence on the wavenumber at different values of the magnetic Rayleigh number.
Here Racm = 6 (black circles), Racm = 10 (open circles) and Racm = 20 (squares).

There are two characteristic lengths in the problem: the magnetic length,lH−2 =and the gravitational length,2(µ∞ − 1)2 (µ2∞ + 1)H∞,
8πηDµ2∞ (µ∞ + 1)

(B 15)1ρg.
Dη

(B 16)

lG−3 =On the scales lH and lG , the magnetic and gravitational Rayleigh numbers,
respectively, are equal to one. Choosing as the characteristic scale of the wavenumber(1ρg/(Dη))1/3 , the dispersion equation in dimensionless form reads (λ̃ = λlG2 /D,
where Ram = (lG /lH )2 is the magnetic Rayleigh number defined by the characteristic
gravitational length, tildes are further omitted)√√√
λ + k2 − k
λ + k2 (2k2 − λ) − 2k3
4 λ + k2− Ram=.
(B 17)kkλ
k2 λ2
In the limit λ → 0, equation (B 17) gives3+ 8k2 .
2k
This gives for the critical magnetic Rayleigh number Racm and wavenumber k∗Ram =Racm = 24k∗2 ,k∗ = (3/32)1/3 .

(B 18)(B 19)The dependence of the growth increment on the dimensionless wavenumber for several
values of the magnetic Rayleigh number is shown in figure 21.
The obtained results give values for the characteristic period of the pattern arisingat the normal field instability and the critical magnetic field strength, which have thesame order of magnitude as observed in the experiment. The critical value of themagnetic Rayleigh number for the critical field in the case of the physical properties

Magnetic micro-convection633

of the magnetic liquid and the density difference of the liquids for the criticalmagnetic field strength and the period of the pattern give the values Hc = 2.47 G
and 2π/k∗ = 27 µm, respectively. These values are close to those observed in theexperiment. As one can see, the formation of the grainy structure on the interface
already appears at field strength of 6.75 G: a value quite close to our estimates. The
period of the structure 20 µm estimated by counting the number of black dots 20
which are seen in figure 14 on the length scale 4 × 10−2 cm also is close to ourestimate.
It should be remarked that scales of this phenomenon are completely different fromthe classical normal field instability, whose scale is determined by the capillary lengthand is of the order of centimetres.
REFERENCESC ABRAL, B. & L EEDOM, L. 1993 Imaging vector fields using line integral convolution. Proc. ACM
93, 263–270.
C EBERS, A. 1981 Dynamics of magnetostatic instabilities. Magnetohydrodynamics 17, 113–121.
C EBERS, A. 1997 Stability of diffusion fronts of magnetic particles in porous media under the actionof external magnetic field. Magnetohydrodynamics 33, 48–55.
C EBERS, A. & M AIOROV, M. M. 1980 Magnetostatic instabilities in plane layers of magnetizable
fluids. Magnetohydrodynamics 16, 21–27.
C HEN, C.-Y. 2003 Numerical simulations of fingering instabilities in miscible magnetic fluids in a
Hele-Shaw cell and the effect of Korteweg stress. Phys. Fluids 15, 1086–1089.
C HEN, C.-Y. & W EN, C.-Y. 2002 Numerical simulations of miscible magnetic flows in a Hele-Shawcell. Radial flows. J. Magnetism Magnetic Mater. 252, 296–298.
D EREC, C., B OLTENHAGEN, P., N EVEU, S. & BACRI, J.-C. 2008 Magnetic instability between
miscible fluids in a Hele-Shaw cell. Magnetohydrodynamics 44, 135–142.
E RGIN, F. G., WATZ, B. B., E RGLIS, K. & C EBERS, A. 2010 Poor-contrast particle image
processing in micro scale mixing. Proceedings of the ASME 10th BiennialConference on
Engineering Systems Design and Analysis, pp. 1–5. ASME.
G OYAL, N. & M EIBURG, E. 2006 Miscible displacements in Hele-Shaw cells: two-dimensional basestates and their linear stability. J. Fluid Mech. 558, 329–355.
I GONIN, M. 2004 Hydrodynamic instabilities of miscible and immiscible magnetic fluids in a
Hele-Shaw cell. PhD thesis, D. Diderot University Paris 7.
I GONIN, M. & C EBERS, A. 2003 Labyrinthine instability of miscible magnetic fluids. Phys. Fluids
15, 1734–1744.
JACKSON, D. P., G OLDSTEIN, R. E. & C EBERS, A. O. 1994 Hydrodynamics of fingering
instabilities in dipolar fluids. Phys. Rev. E 50, 298–307.
M AIOROV, M. & C EBERS, A. 1983 Magnetic microconvection on the diffusion front offerroparticles. Magnetohydrodynamics 19, 376–380.
M EZULIS, A. & B LUMS, E. 2005 Experimental investigations of the microconvective instability in
optically induced gratings. Magnetohydrodynamics 41, 341–348.
ROSENSWEIG, R. E. 1985 Ferrohydrodynamics. Cambridge University Press.
S AMARSKIJ, A. A. & G ULIN, A. V. 1989 Numerical Methods. Nauka.
TAN, C. T. & H OMSY, G. M. 1988 Simulation of nonlinear viscous fingering in miscible
displacement. Phys. Fluids 31, 1330–1338.
W EN, C.-Y., C HEN, C HING -YAO & K UAN, D.-C. 2007 Experimental studies of labyrinthine
instabilities of miscible ferrofluids in a Hele-Shaw cell. Phys. Fluids 19, 084101.
Z IMMERMAN, W. B. & H OMSY, G. M. 1992 Viscous fingering in miscible displacements:
unification of effects of viscosity contrast, anisotropic dispersion, and velocity dependence ofdispersion on nonlinear finger propagation. Phys. Fluids 4 (11), 2348–2359.