Dynamics and stability of three-dimensional ferrofluid films in a magnetic field

We consider the interfacial dynamics of a thin, ferrofluid film flowing down an inclined substrate, under the action of a magnetic field, bounded above by an inviscid gas. The fluid is assumed to be weakly conducting, and its dynamics are governed by a coupled system of the steady Maxwell, Navier–Stokes, and continuity equations. The magnetization of the film is a function of the magnetic field, and is prescribed by a Langevin function. We make use of a long-wave reduction in order to solve for the dynamics of the pressure, velocity, and magnetic fields inside the film. The potential in the gas phase is solved by means of Fourier Transforms. Imposition of appropriate interfacial conditions allows for the construction of an evolution equation for the interfacial shape, via use of the kinematic condition, and the magnetic field. We study the three-dimensional evolution of the film to spanwise perturbations by solving the nonlinear equations numerically. The constant-volume configuration is considered, which corresponds to a slender drop flowing down an incline. A parametric study is then performed to understand the effect of the magnetic field on the stability and structure of the interface.

Free surface flows involving ferrofluids have received some attention in the literature, starting with the linear stability analysis of Cowley and Rosensweig [11] who demonstrated the presence of the so-called normal field instability (NFI). This corresponds to the transition of a magnetic film of initially uniform thickness into an array of hexagonal, crest-like structures in the presence of a magnetic field oriented perpendicularly to the film. The NFI, which occurs beyond a critical value of the magnetic induction, has since been studied by Zelazo and Melcher [12] who showed that normal and tangential magnetic fields are destabilizing and stabilizing, respectively. More recently, Müller [13] showed that viscosity can delay the onset of NFI, and Gollwitzer et al. [14] found good agreement with the work of Cowley and Rosensweig [11] through a Langevin function for the relationship between magnetization and the magnetic field; a similar function will be used in the present paper.
For thin film flows, Joo [15] used the long-wave approximation to derive a nonlinear evolution equation for the ferrofluid film accounting for weak magnetic effects. Moulton and Pelesko [16] used a combination of lubrication theory and experiments to study ferrofluid film drainage. In these models, magnetic effects are accounted for through the inclusion of Maxwell's stresses at the interface and in the bulk. The relation between the magnetic field and the magnetization in the case of ferrofluids is either linear [12] or, for sufficiently large magnetic fields, nonlinear [6]; in the former case, it is possible to rationalise ferrofluid flows by analogy with their electrified counterparts [17][18][19][20][21].
The instabilities exhibited by flows involving thin films and slender drops have been well studied [1,2]. In situations where the film, or drop, flow down an incline, a front develops at the film/drop leading edge, which is unstable to spanwise perturbations, under certain conditions. Investigations of the front dynamics tend to fall into either the 'constant-flux' [22] or 'constant-volume' [23] cases. For the case of a constant-volume of fluid released on an inclined plane, a capillary ridge forms behind a contact line, and eventually becomes unstable to spanwise perturbations. Over time, finger-like structures form, and it was first shown by Huppert [24] that the shape of the patterns can be determined by the contact line dynamics. Also, the presence of the contact line creates a paradox when the no-slip boundary condition is used, leading to multivalued velocity fields [25,26]. To overcome this issue, it is common to apply a slip condition, or include a precursor film ahead of the contact line.
The presence of a ferrofluid in the fingering instability problem will lead to a new effect that has not yet been explored. Due to the magnetic particles that are dispersed within the fluid, a Maxwell pressure is generated following the application of a magnetic field. It is possible to examine the effect of a magnetic field on the fingering instability of a thin film or slender drop using the long-wave approximation [1,2] to couple the interfacial dynamics to the Maxwell stresses. Using this approach, the case of a highly conducting ferrofluid with a nonlinear magnetic susceptibility was investigated by Conroy and Matar [6] who found that the magnetic field has a destabilizing effect and led to a complex surface structure that increased in amplitude with increasing magnetic pressure; this work, however, was restricted to two-dimensional flows.
In this paper, we extend the model developed by Conroy and Matar [6], for the long-wave evolution of a ferrofluid film, to three dimensions. We solve the resultant equations numerically using pseudo-spectral methods, and investigate the spanwise instability of a falling drop down an inclined plane; the magnetic effects enter through a Maxwell pressure. We find that increasing the Maxwell pressure causes the front velocity and the amplitude of the ridge formed at the advancing contact line to increase. We also find magnetic effects to be destablizing leading to finger-formation, with their growth rate increasing with the Maxwell pressure. In addition, the fingers develop spike-like structures that continue to advance downstream from the initial point of deposition. The amplitude of these structures was also found to increase with magnetic field strength.
The rest of the paper is organized as follows. In Sect. 2, we formulate the three-dimensional model for a thin ferrofluid film in a magnetic field, and use long-wave theory to derive a single equation for the interface. In Sect. 3, we provide a discussion of the numerical solutions of this equation, focusing on the stability characteristics of the interface. Finally, in Sect. 4, we provide concluding remarks.

Formulation
We consider the dynamics of a ferrofluid film with magnetization M on a solid, rigid, and impermeable plane, inclined at an angle θ to the horizontal. We use a rectangular coordinate system (x, y, z) to model the film flow in which the wall is located at z = −z w . The interface, which has a constant surface tension γ , is defined by z = S(x, y, t), and separates the film from the surrounding gas. Far away from the film, a magnetic field is applied in a direction normal to the inclined plane. The substrate beneath the film is assumed to be maintained at a constant magnetic field H w . Here, we consider the constant-volume configuration, corresponding to a drop placed on a substrate, with a thin precursor film used to relieve the stress singularity at the moving contact lines.

Governing equations
The governing equations describing ferrofluid films were presented previously in the work of Conroy and Matar [6]; thus, we only provide the dimensionless form, suitably extended to three dimensions. We scale the vertical dimension with the initial film height, z w , as (z, S) = z w (ž,Š), the horizontal dimensions as (x, y) = L(x,y), the fluid velocity as (u, v, w) = V (ǔ,v, δw), the pressure as p =pμV L/z 2 w , and time as t =ť L/V ; here, L is a characteristic lateral extent, V is a scale for the velocity to be defined later, μ is the fluid viscosity, and δ ≡ z w /L 1 is a small parameter, which will be used as the basis for the asymptotic reduction to be carried out below. The magnetic film in each layer is scaled as H i = H ∞Ȟi , the magnetic potential as φ i = H ∞ z wφi , and the magnetization as M i = H ∞Mi , where i = 1 and 2 for the film and gas, respectively, and H ∞ is the incident magnetic field. Henceforth, we drop the check decoration for the dimensionless quantities.
In the film (contained in a domain defined by −1 ≤ z ≤ S) and surrounding gas (z ≤ S), the magnetic field can be defined in terms of the potential as H i = −(δφ i x , δφ iy , φ iz ), where the magnetic potential, φ i , is determined from in the film and surrounding gas, respectively. The parameter β F is related to the dimensionless magnetization defined from the relationship where β = M s /H ∞ , ξ = ξ H ∞ = 3χ 0 H ∞ /M s , χ 0 is the initial susceptibility, M s is the saturation magnetization, For the free surface, z = S(x, y), we define the normal and tangential vectors as well as the normal and tangential jump conditions from Eqs. (1) and (2) in dimensionless form: At the wall, z = −1, we fix the magnetic field as was done by [27], which is different from the model of [6]. The appropriate condition in this case for a flat surface is where we have set the substrate magnetic field to a small value so that it appears in the leading-order terms to be discussed in the next section. Far from the film, the magnetic field is H 2 → H ∞ as z → +∞.
The dimensionless momentum and continuity equations are Re Re where , and the dimensionless Maxwell pressure reads Here, the dimensionless groups appearing in the momentum equation are defined as and represent the Reynolds number wherein ρ is the fluid density, a magnetic parameter in which μ 0 (not to be confused with the viscosity) is the magnetic permeability, and a gravitational parameter, respectively. At the interface, the normal and tangential stress balances are where Ca = μV /γ is the capillary number, γ is the surface tension, and K is the curvature of the interface: and the fluid stress is Finally, the dimensionless kinematic condition is and, at the substrate surface, z = −1, u = v = w = 0 corresponding to the no-slip and no-penetration conditions.

Long-wave approximation
We seek a solution to the magnetic potential for a large effective magnetic permeability of the ferrofluid. In this way, the film, as viewed from the outer field, resembles a magnetized sheet. With this assumption, we take the order of magnitude of the dimensionless groups as β = O(1) and ξ = O(1) − O(10) as discussed in [6]. The surrounding gas region is not slender so we re-scale as z =zδ −1 [19]. We now set the velocity scale to V = δ 3 γ /μ so that the capillary number in Eq. (16) is Ca = μV /γ = δ 3 , and the capillary pressure term in this equation is of order unity. We also setQ H = δ 2 Q H so the Maxwell pressure is of the same order as the capillary pressure. In addition, we take ξ =ξ/δ, F =Fδ −1 , and H ∞ = 1 + δ H 1 (x). Therefore, the magnetization term is expanded as and Eqs. (1) and (2) for the magnetic potential in the ferrofluid film and gas, respectively, become The far-field boundary condition is and the interfacial conditions corresponding to Eqs. (5), (6) are expressed as The substrate conditions corresponding to Eqs. (7), (8) are now re-written as where the solid wall is at z = −1.
We now expand the potential in a perturbation series of the form φ i = φ i0 + δ 2 φ i1 + · · · (for i = 1 and 2), which is asymptotically valid in the limit of δ → 0, and collect terms of equal powers of δ. The idea is that we are determining the linear solution and its perturbation in powers of a small parameter δ 2 (for more details, see Hinch [28]). From Eq. (20), the potential to leading-order (i.e. the equation with order δ 0 ) is βFφ 10z = c(x, y) in the film, where c(x, y) = 0 from the boundary conditions. Integrating again, the potential is only a function of x and y, φ 10 = φ 10 (x, y), and is determined from the next order of the expansion. Also to leading-order φ 10 which is only a function of x and y.
For the order δ 2 equation, the potential in the film to the next order is determined from the following set of equations: Integrating Eq. (28) once with respect to z, inserting into the above boundary conditions, and subtracting we get which is a condition for the horizontal gradient of the potential. The potential in the outer gas is more easily determined through subtracting the outer magnetic field by defining with the following far-field and interfacial conditions, respectively, given by We solve Laplace's equation in the gas phase by using the Fourier transform defined aŝ wherek x andk y are the wavenumbers in the x and y-directions, respectively. Here, the integrations are taken to infinity in x and y; however, it is possible to integrate the function Φ over a finite domain provided it is periodic at the boundaries. In this analysis, we assume that the surface height, velocity, pressure and magnetic potential have periodic boundary conditions in the x-y plane, and perform the integrations over a finite domain. Laplace's equation in Fourier space becomeŝ wherek 2 =k 2 x +k 2 y . Using the boundary conditions far from the film, Φ 2z = 0 asz → +∞, and the condition given by Eq. (34), the potential in spectral space iŝ where S H ∞ is the Fourier transform of the product of S and H ∞ and the derivative of Eq. (38) at δS = 0 iŝ In the dimensionless momentum equations, we take Re = O(1) and expand in powers of δ (e.g. u = u 0 + δu 1 + · · · ). The continuity and momentum equations to leading order are whereĜ = δG cot(θ ) and At the interface z = S(x, y), the normal and tangential jump conditions to leading-order are where (H 20 · n − H 10 · n) 2 = δ −1 + 2(H 1 (x) + Φ 20z ) and δ have been retained here for completeness; however, it will be eliminated later when we take the derivative of this quantity. Integrating the leading-order momentum equations, and using the interfacial (air-film interface located at z = S(x, y, t)) and no-slip boundary conditions at z = −1, yields the following expressions for the film pressure, and streamwise and spanwise velocity components, respectively: From the kinematic condition, the interface evolves according to where the two-dimensional function K (x, y) is Here, we have defined H ∞ = 1 + δ H 1 (x). Note that the presence of the O(δ −1 ) in the expression for K (x, y) is inconsequential since K x and K y (rather than K (x, y)) are required to form u and v, respectively. In Eq. (52), the vertical derivative of the magnetic potential is where the potentialφ 10 comes from Eq. (31) transformed to Fourier space.

Results
The governing equations describing the flow of a thin ferrofluid in a magnetic field are given by Eqs. (31), (39), (51), and (52). The dimensionless groups governing the process are: the saturation magnetization relative to the outer magnetic field, β, the magnetic parameter,Q H , the scaled initial susceptibility,ξ , and the gravitational numbers G, andĜ = δG cot(θ ). For typical laboratory conditions, using magnetite as the suspended particles, we can assume the following parameter values: M s ≈ 1.5 × 10 4 A/m, ξ ≈ 2 × 10 −4 − 5 × 10 −4 m/A, and H ∞ ∼ M s [14]. The order of magnitude for the dimensionless groups are then: The effect of the magnetic and fluid dynamics parameters were studied by [6]. Here, we focus on the magnetic field strengthQ H , and fix the other parameters to β = 1,ξ = 5, G = 0.50, andĜ = 0.25.

Linear stability analysis
It is instructive to analyse the linear stability for a flat interface in order to gain insight into the dynamics of a ferrofluid film near instability onset. We do this by linearizing the governing equations about the uniform base state S = 0 and φ = φ b (a constant) using normal modes of the form: Here, k x and k y are the wavenumbers in the x and y-directions, respectively;S andφ are constants to be determined below; and S and φ are the perturbations about the base state. Introducing the expansions into Eqs. (31) and (39) and retaining terms on the order of the perturbation, the magnetic potential is expressed as with the magnetic potential being linearly related to the interface height and k 2 = k 2 x + k 2 y . We take H w = 0 for simplicity, and from Eqs. (51) and (52), the real part for the dispersion relation is obtained as 3ω r =Q H k 3 The first term is positive so the magnetic field is destabilizing, whereas capillary forces always stabilize the film for sufficiently large k as found by [6,11,29]. The dispersion relation, for k y = 0 and k x = k, is displayed in Fig. 1, showing the wavenumber band for instability. The result is similar to the electrohydrodynamic case [17], where the electric field, for a perfectly conducting viscous film, was also found to be destabilizing. In the large βξ k limit, we simplify the dispersion relation as ω =Q H k 3 /3 −Q H k 2 /βξ − k 4 /3 −Ḡk 2 /3. Taking the derivative of this equation with respect to k and setting it to zero, we obtain the location of the maximum growth rate as In Fig. 2, the maximum growth rate and the most dangerous mode are shown as a function ofQ H and β. The most dangerous mode exhibits a jump followed by a linear increase withQ H . The large βξ k solution from Eq. (57) is in close agreement with the full solution even for moderate values of β. The prediction of the jump in k max from Eq. (58) is also very good for sufficiently large βξ k.

Numerical results
We solve the governing equations using pseudo-spectral methods with a sufficiently large number of modes in the x-y direction. Typically, we use 256 × 64 modes in the simulations without any discernible difference in the results with half or double the number of modes. Two-dimensional Fast Fourier Transforms are used for the basis functions and the boundary conditions are periodic. The thin film equations are discretized in time using an explicit RK-3 method, and the magnetic field equations are solved using an iterative method. Domain decomposition is used to parallelize the code using MPI-2. With this technique, the spatial domain is broken up into strips along the y-direction with each processor computing n y /n p × n x points, where n y is the number of modes in the y-direction, n p the number of processors and n x the number of modes in the x-direction. The main MPI communications occur in the FFT routines in order to perform the two-dimensional integrations. We checked that mass conservation is maintained in the computations and that the solution remains unchanged with the increasing resolution. A precursor film is used to avoid the contact line stress singularity with a thickness that is sufficiently small. For an initial profile, we use a Gaussian distribution: Here, S f = −0.98 (the precursor film thickness is then 0.02) is the location of the precursor film, and we take x 0 = 5 unless otherwise specified. Given the initial surface profile, we obtain the initial magnetic potential by solving Eq. (31). In addition, we use periodic boundary conditions at the ends of the x-y domain for S(x, y, t) and φ 10 (x, y, t).

Two-dimensional unperturbed case
First, we investigate the case of a two-dimensional (2D) flow in the absence of imposed spanwise perturbations, which will be used in the subsequent stability analysis to be discussed below. In Fig. 3, we show snapshots of the spatial development of the interface S at t = 80 for three different values ofQ H , including the non-magnetic case, withξ = 5, β = 1, G = 0.5, andĜ = 0.25 held constant. For all values ofQ H considered (includingQ H = 0, Fig. 4 Parameteric dependence of the steady-state value of S on Q H for a three-dimensional axisymmetric drop spreading on a flat substrate in the presence of a magnetic field. Here, ξ = 5, β = 1, and G = 0.0 =Ĝ = 0.0. We also use the same initial Gaussian profile for the surface height with x 0 = π as expected), the drop forms a capillary ridge near the advancing 'contact line' located at the downstream leading edge, and a long tail which extends back to the initial drop location. The magnetic field creates a jump in φ 10x , causing an increased pressure at the peak. As a result, the drop accelerates and the capillary ridge steepens and increases in height. In addition, for the largest values ofQ H examined, a 'secondary' ridge appears in the tail region with a trough forming between this and the 'primary' ridge at the leading edge. To isolate the effect ofQ H (from those associated with the gravitational acceleration) on the primary ridge structure, we examine briefly the case of a drop on a flat substrate in the presence of a magnetic field. The drop structure is initialized by the same Gaussian profile as that used to generate the results in Fig. 3, and the simulation is carried out until a steady-state solution is reached. As shown in Fig. 4, for increasing values ofQ H , the peak height of the drop, S max , increases, and the drop narrows due to mass conservation, forming spike-like structures.
The time evolution of the maximum height of the drop, S max , and that of the spatial location corresponding to S max , x h , are plotted in Fig. 5. For the non-magnetic case,Q H = 0, S max decreases with time starting from t = 0 as it spreads down the inclined plane. When the magnetic pressure is applied, i.e. forQ H > 0, S max increases initially as the surface is pulled in the direction of the incident magnetic field but later decreases as the constant volume drop flows down the incline. It is also evident from Fig. 5 that the spreading is accelerated via the application of the magnetic field through inspection of the temporal evolution of x h which increase withQ H for all t values considered.
A semilogarithmic plot of x h vs. t shown in Fig. 3c reveals that the spreading of the drop exhibits a power-law behaviour with a well-defined exponent at large times. In all of the cases shown, the exponent is similar for different values ofQ H although the slope is a little steeper for larger vales of the magnetic pressure. It appears from the power-law exponent that the motion of the ridge down the incline is dominated by gravity despite the presence of the magnetic field for all theQ H investigated. Neglecting magnetic and capillary effects in Eq. (51), and taking into consideration spreading in the x−direction only, we have S t + G S + 1) 3 /3 ∼ 0. The continuity constraint then leads to x h ∼ t 1/3 , which is in agreement with the numerically determined exponents (see Fig. 5).

Single spanwise mode
We now examine the stability of the 2-D states discussed in the foregoing. We perform the stability analysis by perturbing a 2-D simulation at time t p , which occurs at a time where the front is sufficiently developed. The front position is perturbed by applying a spanwise disturbance, which is added to the right-hand-side of Eq. (51) expressed by where A 0 is the amplitude, k p is the wavenumber, and G a (t) is a Gaussian function of time that is centred at t = t p . In the absence of perturbations, the front is stable for very long times, as has been verified numerically. The simulations are performed in a moving reference at the speed v f , which was shown by [22] not to affect the results. We also note that in this case, the front is stable to spanwise perturbations in the absence of a Maxwell pressure. In order to avoid the complications associated with the temporally evolving base state in the linearized form of the governing equations, we determine the rate of growth of the spanwise stability numerically. To do this, we define a scaled spanwise 'energy', E t , and a corresponding growth rate, ω t , as This expression is similar to that used by [30], except that we are using the derivative of the interface in the spanwise direction S y to represent the growth in energy of the perturbation. It is important to note that our approach does not correspond to a linear stability analysis: we have not linearized Eq. (51), and indeed the amplitude of the imposed perturbation given by Eq. (59), although initially small, is permitted to grow such that the disturbance can enter into the nonlinear regime. Thus, ω t does not correspond to the growth rates that can be obtained from a linear stability analysis (except perhaps over a limited time duration in which the perturbation has grown exponentially in amplitude yet its amplitude has remained sufficiently small so that we are still in the linear regime). The values of E t and ω t are tracked during the course of the simulation from t p onwards for a given value of k p , as shown in Fig. 6a, b. It is seen clearly that E t and ω t undergo an increase with time for t > t p indicating disturbance growth, with both growth measures increasing withQ H , which highlights the destabilizing effect of the imposed magnetic field. From Fig. 6b, it is evident that following the initial rapid increase, the disturbance growth rate decreases at later times; in the case ofQ H = 2.5, ω t appears to asymptote to a constant value.
In Fig. 6c, we plot numerically constructed 'dispersion relations' showing the dependence of ω t on k p at t = 30 forQ H = 2.5, 3, with all other parameters remaining unchanged from those used to generate Fig. 3. It is seen clearly that the dispersion curves are structurally similar to those shown in Fig. 1, associated with the linear stability of a film of uniform thickness, and exhibit low-and high-wavenumber cut-off 'modes', as well as a 'most dangerous mode' at an intermediate k p value. Increasing the value ofQ H leads to an increase in the wavenumbers associated with both the cut-off and most dangerous modes, as well as a rise in growth rate for the range of 'unstable' wavenumbers. This further highlights the fact that magnetic effects exert a destabilizing influence on the dynamics. In Fig. 6d, we also show the effect of the perturbation start time, t p , on ω t forQ H = 3.0, and t = 25, 30, 35. For the range of t p values shown, ω t decreases with t (see also Fig. 6b) but appears to be relatively insensitive to t p for t p > 10. This result shows that the growth rates are not that sensitive to the perturbation start time even though the two-dimensional base state is evolving as the constant volume fluid spreads out and thins. Because the unperturbed flow is evolving with time there is a strong dependence in ω t with time as shown in Fig. 6b. This result is in contrast with the constant flux configuration, where the unperturbed film reaches a steady state and the initial growth rate is a weak function of time.
In Fig. 7, we plot the spatiotemporal evolution of the interface associated with the application of a spanwise perturbation. At t = t p = 10, the drop has already spread down the incline, forming a capillary ridge at the downstream leading edge, and has just interacted with a spanwise, cosine perturbation of k p = 1 in its path. Due to the Maxwell pressure, the perturbation grows, as indicated by the results presented in Fig. 6, giving rise to a pronounced hump at the leading edge of the spreading drop. Since the magnetic field is oriented normal to the substrate, the peak of the hump is pulled in the wall-normal direction and a spike-like structure forms (see the surface plot associated with t = 100 in Fig. 7). With the increasing time, the perturbed interface acquires the shape of a long finger, and the hump decreases in amplitude due to mass conservation. The large decrease in the amplitude  of the finger over time is due to the spreading of a constant volume of fluid. In the constant-flux configuration, we expect to get a similar finger structure except that the height will probably not decrease as much over time.
Surface plots (top views) of the interface for different spanwise perturbation wavenumbers are shown in Fig. 8. Here, k p = 0.5, 1, 1.5,Q H = 3, and the rest of the parameters remain unchanged from Fig. 3. Since the domain size was chosen to fit one spatial period, only one complete finger is shown in each of the panels. The diameter of the head appears to have a similar size, in the x-y plane for all three k p values. In addition, the shoulder of the front is smaller for increasing values of k p . The growth of the finger length with time, defined as the streamwise distance between the thickened finger hump and the finger root, is shown in Fig. 9 for k p = 0.5, 1 and 1.5 andQ H = 3. The fastest growing mode at early times occurs for k p = 1.5, which is consistent with the results presented in Fig. 6c. The finger growth at larger times is clearly much larger for k p = 1, and we see an acceleration phase for intermediate values of time, followed by a slower period that is likely due to the decrease in volume in the hump of the thickened finger.

Multiple spanwise modes
Next, we investigate the stability of the film to a large set of spanwise modes. In this way, the most dangerous mode will be naturally selected with the higher wavenumber content dampened by capillarity. In reality, the distribution of the perturbation will be non-uniform and environmentally dependent, but since the distribution is unknown, and we are interested in the self selection process, we use a uniform distribution. Also, we initiate the instabilities after a time t p following [22], since the time at which noise is imparted on the fluid is unknown.
Surface plots of the interface are shown at t = 150 in Fig. 10 forQ H = 2.5, 3, 3.5. For small values ofQ H , the interface profile looks similar to the non-magnetic case, consisting of a relatively smooth profile with thickened ridges at the leading edge. For large values ofQ H , the drop velocity increases, and the tips of the fingers develop into spike-like structures. Also, for increasingQ H , the number of fingers increase indicating that the most dangerous wavelength decreases with a stronger magnetic field.

Conclusions
In this paper, we have investigated the dynamics of a thin drop of a ferrofluid spreading down an inclined plane in the presence of a magnetic field oriented in a wall-normal direction. We have developed a model of the flow using a long-wave analysis for a viscous fluid accounting for Maxwell stresses. The magnetic field equations were solved for a relatively large magnetic permeability, and a magnetization that was assumed to be a nonlinear function of the magnetic field using a spectral technique. The stress singularity at the moving contact lines was relieved through the use of a thin precursor layer. For two-dimensional flow, we have shown that the spreading process is accompanied by the formation of a thickened ridge at the advancing contact line. With the increasing magnetic field strength, another hump is seen to develop upstream of the ridge. The spreading is observed to follow a power-law, which is consistent with gravity-driven spreading for the range of parameters investigated in the present study.
We have also examined the stability of the interface to the applied perturbations. We have derived a dispersion relation for the case of a flat ferrofluid film, which demonstrates the destabilizing effects of the magnetic field in the linear regime. We have also studied the stability of the two-dimensional spreading drops to spanwise perturbations of either single or multimodes, and shown that significant growth takes place leading to the formation of long fingers with spike-like ridges; the latter become particularly pronounced with increasing magnetic field strength. The results of this work will shed light on how magnetic fields can be used to manipulate fluids, which will be relevant to microfluidics and lab-on-a-chip applications.