# Simulation of Thermomagnetic Convection in a Cavity Using the Lattice Boltzmann Model

Authors: Mahshid Hadavand; Antonio C. M. Sousa
Journal: Journal of Applied Mathematics
 (2011)

---

## Abstract

Thermomagnetic convection in a differentially heated square cavity with an infinitely long third dimension is numerically simulated using the single relaxation time lattice Boltzmann method (LBM). This problem is of considerable interest when dealing with cooling of microelectronic devices, in situations where natural convection does not meet the cooling requirements, and forced convection is not viable due to the difficulties associated with pumping a ferrofluid. Therefore, circulation is achieved by imposing a magnetic field, which is created and controlled by placing a dipole at the bottom of the enclosure. The magnitude of the magnetic force is controlled by changing the electrical current through the dipole. In this study, the effects of combined natural convection and magnetic convection, which is commonly known as “thermomagnetic convection,” are analysed in terms of the flow modes and heat transfer characteristics of a magnetic fluid.


---

## Body

## 1. Introduction

Natural convection in enclosures is of interest in many engineering applications, such as the cooling of electronic equipment and solar energy collection, and it has been extensively studied experimentally and numerically by several researchers [1–4]. The enhancement of convective heat transfer in these devices is an active area of research. For small length-scale applications, such as the cooling of microelectronic devices, thermomagnetic convection has been identified as a viable approach for augmenting and controlling the convective heat transfer. This involves controlling the motion of colloidal suspensions containing magnetic nanoparticles, known as ferrofluids, using external magnetic fields [5, 6].The study of flow and heat transfer characteristics of natural convection of a magnetic fluid has received considerable attention lately. Finlayson [7] in 1970 explained how an external magnetic field imposed on a ferrofluid with varying magnetic susceptibility, for example, due to a temperature gradient, results in a nonuniform magnetic body force, which leads to thermomagnetic convection. This phenomenon is complex in nature, and so a thorough understanding of the relation between an applied magnetic field and the resulting heat transfer is necessary for the proper design and control of thermomagnetic devices [8].The numerical studies of heat transfer with magnetic fluids reported in the literature are relatively scarce and invariably use as governing equations the Navier-Stokes and energy equations. In what follows, relevant work in the field is succinctly reviewed. Yamaguchi et al. [9] experimentally and numerically studied the heat transfer characteristics and flow behavior for natural convection of a magnetic fluid in a rectangular cavity under an imposed uniform vertical magnetic field. Al-Najem et al. [10] numerically calculated the flow and temperature fields under uniform magnetic field in a tilted square cavity with isothermal vertical walls and adiabatic horizontal walls at a low Prandtl number of 0.71 and demonstrated that the suppression effect of the magnetic field on convection currents and heat transfer is more significant for low inclination angles and high Rayleigh numbers. Krakov and Nikiforov [11] determined the influence of the angle between the direction of the temperature gradient and that of a uniform magnetic field on the convection structure and intensity of heat flux in a square cavity.A few other researchers have considered spatially nonuniform magnetic fields for the experimental or numerical investigations but did not fully describe the variations in those fields [12–15]. For example, Sawada et al. [12] and Kikura et al. [13] experimentally investigated the influence of a varying magnetic field in a cubical enclosure and in concentric horizontal annuli, respectively, but the field gradient was only measured in a single direction. Moreover, the orientation of the magnetic field was not described. Hadavand et al. [16] simulated thermomagnetic convection in a cavity which was in close proximity to a line dipole. The imposed nonuniform magnetic field was considered only in one direction and the effect of temperature on the magnetic susceptibility was ignored. Ganguly et al. [8] performed a study to simulate thermomagnetic convection in a square enclosure by introducing a line dipole using a Navier-Stokes solver. The effect of temperature on the magnetic susceptibility was considered; because they used a very strong nonuniform magnetic field, the buoyancy was neglected in their simulation. They found that the average Nusselt number on the wall increases with increasing the magnetic dipole strength and temperature but decreases by increasing the fluid viscosity.The magnetic control has enabled numerous developments dealing with technical and medical applications, and these suspensions are already being considered the next generation heat transfer fluids as they offer the possibility of achieving heat transfer rates much higher than those of conventional heat transfer fluids and fluids containing microsized metallic particles [17].In the present work, the objective is to numerically study the effects of combined natural convection and magnetic convection in a square enclosure using the lattice Boltzmann method. The square enclosure is under the influence of an imposed nonuniform two-dimensional magnetic field that satisfies Maxwell’s equations. The variation of the ferrofluid susceptibility with the temperature difference is considered in this study.


## 2. Physical Model

In the absence of an external magnetic field, the nanoparticles are oriented randomly in the ferrofluid; once the external magnetic field is applied, the nanoparticles align with the applied magnetic field. The generated force between the magnetic field and the homogenously distributed magnetic nanoparticles enables the manipulation of the ferrofluid by adjusting the applied external magnetic field. This uniformly distributed magnetic force manifests itself as a body force analogous to the gravity field. These suspensions, in general, exhibit normal liquid behavior coupled with super paramagnetic properties. This leads to the possibility of controlling the properties and the flow of these liquids with relatively moderate magnetic field strengths. Each solid magnetic nanoparticle in the ferrofluid is covered with a dispersant material in order to keep the particles separated, and, in this way, preventing agglomeration of the particles due to gravity, coagulation, or an external magnetic field.Figure1 presents a schematic of the enclosure used for the numerical simulation, which has height D and its vertical walls are at constant, but different, temperatures. The third dimension is considered to be infinite; therefore, the flow in this configuration can be assumed to be two-dimensional. The left and right vertical walls are at constant high temperature, Th, and constant low temperature, Tc, respectively. The upper and lower walls are adiabatic. A line dipole, which provides the external magnetic field, is placed adjacent to the lower adiabatic wall at a distance d. The magnetic field is governed by the Maxwell relations in static form [18], as follows:(2.1)∇·B=0,∇×H=0,
where B is the magnetic field inside the ferrofluid due to the line dipole which can be expressed as follows [8]:(2.2)B=μ0(1+χ)m[sinθr2er-sinθr2eθ],
where H and B are related through the following equation:(2.3)H=1μ0(1+χ)B,
where μ0 is the magnetic permeability of free space, μ0=4π×10-7(N/A2); m denotes the magnetic dipole moment of the electromagnet coil (having two parallel counter-conductors separated by a distance b and carrying current I) per unit length (m=Ib/2π). The magnetic fluid susceptibility, χ, varies with temperature as follows:(2.4)χ=χ01+β(T-T0),
where T0 is the reference temperature and χ0 is the magnetic field susceptibility at reference temperature.Figure 1
Schematic of the cavity and position of the line dipole magnet.The magnetization of the ferrofluid occurs when it becomes polarized in the presence of an external magnetic field. Consequently, attractive forces acting on each particle are produced due to the interaction of the ferrofluid polarization and the external magnetic field. The attractive force acting on the particles can then be treated as a body force acting on the ferrofluid, which is dealt with as a homogenous medium. The attractive force on the ferrofluid per unit volume is(2.5)F=(M·∇)B,
where M is the magnetization vector and is related to H by M=χH. The magnetic field strength diminishes with increasing distance from the dipole. The Kelvin body force, see (2.5), can be written as follows:(2.6)F=μ0χ(1+χ)(H·∇)H=12μ0χ(1+χ)∇(H·H)+μ0χH(H·∇χ).The first term on the right-hand side of (2.6) has a same nature as the pressure term in the momentum equation. The second term becomes important when there is a spatial gradient in the fluid susceptibility. In the absence of a temperature gradient, this term is zero. The Kelvin body force is simplified in another study [19] to the following format by defining an effective pressure as (2.7)P*=P-μ0χ02H2,
and it yields (2.8)F=12μ0χ0[1-β(T-T0)]∇(H·H)+μ0χ02β(H·∇T)H.In this work, the single relaxation time lattice Boltzmann method is used to simulate the thermomagnetic convection in a two-dimensional square enclosure.


## 3. Numerical Model

The lattice Boltzmann method (LBM) is a discrete particle-based mesoscopic approach for viscous fluid flow simulation. The LBM has been shown to be an efficient tool for flow simulation in complex geometries, when compared with conventional fluid dynamics approaches [20]. In this method, finite volumes of fluid are presented as discrete particles and the fluid dynamics is obtained based on a probability distribution function of the particles and a discrete set of prescribed velocities.The thermal lattice Boltzmann model, which is employed in the present work, is a double-distribution-function thermal model. Two distinct distribution functions for the flow and temperature fields aref and g, respectively, where f is the density (or fluid particles) distribution function and g is the temperature distribution function. This present thermal lattice Boltzmann model uses the BGK approximation [21] and it is based on the work of He et al. [22] and Shi et al. [23]:(3.1)fi(x+ciδt,t+δt)-fi(x,t)=-δtτν[fi(x,t)-fieq(x,t)]+fg+F,gi(x+ciδt,t+δt)-gi(x,t)=-δtτc[gi(x,t)-gieq(x,t)],
where δx and δt are the lattice grid spacing and the lattice time step, respectively. These two are related to each other through the streaming speed, c, as c=δx/δt. However, in most of the LBM simulations on a uniform lattice, including those of the present work, δx and δt for simplicity are considered to be equal to one, and fg is the buoyant body force term, which is formulated by using the Boussinesq approximation [24], namely:(3.2)fg=3ωiβρgy(T-Tave)ciy,
where β is the thermal expansion coefficient, gy is the acceleration of gravity acting in the y-direction of the lattice links, and ciy is the y-component of ci. The terms ρ and T are the local density and temperature and are calculated at each lattice site. Tave is the average temperature of hot and cold wall. The Kelvin body force, F, is the force that a magnetic fluid experiences in nonuniform magnetic field that can be calculated based on (2.5).The most common two-dimensional lattice is the D2Q9 lattice, which is used in this study, where “2” denotes the number of space dimensions and “9” refers to the discrete set of nine velocities as proposed by Qian et al. [25]. The D2Q9 lattice and its labelling for the lattice directions are shown in Figure 2. For the D2Q9 lattice, the discrete velocity set is formulated as follows:(3.3)ci=c[cos((i-1)π2),sin((i-1)π2)],i=1,2,3,4,ci=2c[cos((2i-1)π4),sin((2i-1)π4)],i=5,6,7,8,c9=0.Figure 2
Topology of the D2Q9 lattice and labelling of the directions.The corresponding equilibrium distribution functions for the density and temperature distribution functions arefieq and gieq, respectively. These two parameters are in the form of a second-order truncated expansion of the Maxwell-Boltzmann equilibrium function:(3.4)fieq=ρωi[1+1cs2ci·u+12cs4(ci·u)2-12cs2u·u],gieq=ωiT[1+1cs2ci·u+12cs4(ci·u)2-12cs2u·u],
where cs is the lattice sound speed and is equal to 1/3 for a two-dimensional lattice with 9 velocities (D2Q9). ωi is the weighting factor and is defined as follows for the D2Q9 lattice:(3.5)ωi={49,i=9,19,i=1,2,3,4,136,i=5,6,7,8.The macroscopic density(ρ), pressure (P), velocity (u), and temperature (T) on each lattice site are calculated as follows:(3.6)ρ=∑ifi,P=ρcs2,ρui=∑ici,T=∑igi.The kinematic viscosity and the thermal diffusivity are given by(3.7)ν=(τν-12)cs2,α=(τc-12)cs2.The no-slip solid-wall boundary condition is set for the vertical and horizontal walls of the enclosure by using the bounce back method. As mentioned before, the top and bottom walls are adiabatic, and the vertical walls are isothermal, with the left wall at higher temperature than that of the right wall.


## 4. Results and Discussion

The physical domain considered in this study is a two-dimensional square cavity. The parametersD and d are depicted in Figure 1, where D is set to 2 mm and d is equal to D/2. The density and dynamic viscosity of the ferrofluid at 300 K were set equal to 1180 kg/m3 and 1.0×10-3 kg/ms, respectively. The Prandtl number, Pr, and fluid compressibility, β, are taken as 5.5 and 5.6×10-4K-1, respectively.A uniform lattice system is employed for the simulations. The numerical solution should be independent of the lattice size, and to this purpose, three different lattice sizes of81×81, 101×101, and 121×121 were examined; results with a lattice density of 101×101 lattice units were found to be adequate to establish lattice-independent solutions for the range of the parameters used in the present study. In order to validate the thermal LBM predictions, laminar natural convection of air, Pr=0.71, in a cavity was simulated for Rayleigh numbers of 104 and 105. The predicted natural convection results with the 101×101 grid system for Ra=104 and 105, respectively, were in excellent agreement with the benchmark steady flow patterns and temperature fields reported by de Vahl Davis [26]. In Table 1, the predicted average Nusselt numbers on the hot wall are compared against the values reported by de Vahl Davis [26] for the same conditions and, as it can be observed, are in excellent agreement.Table 1
Comparison of the predicted average Nusselt number for two different Rayleigh numbers against the values reported in [26].
Rayleigh numberAverage Nusselt number104Present work2.247de Vahl Davis [26]2.243105Present work4.481de Vahl Davis [26]4.519The average Nusselt number on the side walls is determined from the following relation:(4.1)Nu¯=-1ΔT∑1DdTdx.The two dimensionless parameters which characterize the natural convective flow, Rayleigh number (Ra) and the Prandtl number (Pr), are defined as follows:(4.2)Ra=gyβΔTD3να,Pr=να,
where gy, β, ΔT, and D are the gravitational acceleration, working fluid thermal expansion coefficient, temperature difference between the two vertical walls, and the characteristic length (equal to the width of the cavity), respectively.Presenting the magnetic dipole will disturb the flow motion in the cavity; consequently, the temperature distribution will be affected along with the heat transfer rate at the walls of the cavity. The isotherms for three different values ofmequal to 0.0, 0.05, and 0.1 Am, are shown in Figures 3(a)–3(c) for ΔT=30K (Ra=104). The influence of the induced magnetic field on the temperature distributions is apparent. If the magnetic field magnitude is relatively weak, the buoyancy effect is dominant and the isotherms and streamlines are similar to those of pure natural convection in the absence of a magnetic field. For a stronger magnetic field the isotherm stratification in the core will change, as shown in Figure 3(c); this is due to the suppression of the free convective flow by the magnetic field. By increasing the magnetic field strength, the effect of the natural convection is eventually eliminated and the flow motion of the ferrofluid is governed by the strength of the magnetic field.Isotherms for the thermomagnetic flow in a square cavity forΔT=30K, Ra=104: (a) m=0.0, (b) m=0.05, and (c) m=0.1Am.
(a)(b)(c)Figures4(a)–4(c) report the isotherms at ΔT=60K (Ra=2×104) for m=0.0, 0.05, and 0.1 Am. From these figures, it can be observed that for large values of ΔT, the effect of increasing the magnetic strength on the temperature distribution in the enclosure is more prominent compared to that in the case of a smaller ΔT (Ra=104, Figure 3). This is because the magnetic susceptibility depends on temperature and influences the motion of the ferrofluid in the cavity. In these circumstances, the colder fluid close to the lower temperature wall moves down towards the region of larger field strength due to its higher magnetic susceptibility, while the warmer fluid with lower susceptibility is displaced away from the line dipole and is pushed upward along the hot wall.Isotherms for the thermomagnetic flow in a square cavity forΔT=60K, Ra=2×104: (a) m=0.0, (b) m=0.05, and (c) m=0.1Am.
(a)(b)(c)For sufficiently large values ofm, the natural convection is suppressed and the magnetic field effect is the dominant factor. As a result, the isotherms will be almost independent of Ra number for a relatively high magnetic field strength. The streamlines for the thermomagnetic flow in the square cavity for ΔT=30K and ΔT=60K and for different values of m, 0.0, 0.05, and 0.1 Am, are presented in Figures 5(a)–5(c) and 6(a)–6(c), respectively.Streamlines for the thermomagnetic flow in a square cavity forΔT=30K, Ra=104: (a) m=0.0, (b) m=0.05, and (c) m=0.1Am.
(a)(b)(c)Streamlines for the thermomagnetic flow in a square cavityΔT=60K, Ra=2×104: (a) m=0.0, (b) m=0.05, and (c) m=0.1Am.
(a)(b)(c)Figure7 shows the predicted values of the average Nusselt number on the hot wall for a range of temperature difference, 15≤ΔT≤90, 5×103≤Ra≤3×104, and 0≤m≤0.2. It can be observed that the average Nusselt number increases with both magnetic field strength and temperature difference. The average Nusselt number increases more than double for m=0.2Am compared to the Nusselt number for free convection only.Figure 7
Average Nusslet number for different temperature difference and magnetic field strength.To show the effect of magnetic field strength on heat transfer enhancement, the buoyancy force was turned off and the average Nusselt numbers on the hot wall were compared with those from the combination of buoyancy and magnetic force. In Figure8, the average Nusselt numbers are presented for different magnetic field strengths and a range of temperature difference. For magnetic field moment, m, lower than 0.1 Am the buoyancy force has a noticeable effect on the Nusselt number, but with increasing magnetic field strength the effect of the buoyancy force is reduced and the Nusselt number can be calculated by considering only the magnetic force for m larger than 0.15 Am.Figure 8
Average Nusslet number based on the magnetic force and combination of magnetic and buoyancy force. (The solid lines are the Nuave based on the combination of the magnetic and buoyancy force. The marked lines are calculated based on the magnetic force only.)


## 5. Conclusion

Simulations for thermomagnetic convection were conducted for a range of Rayleigh numbers and magnetic field strengths in the laminar regime. The effects of these parameters on the heat transfer characteristics were analyzed; the numerical results showed that with the application of an external magnetic field, the temperature and velocity fields were significantly modified.For weak magnetic fields, both buoyancy force and magnetic force have effects on the heat transfer. However, for sufficiently high magnetic field strength, the magnetic convection is suppressed for all the examined values of the Rayleigh number. This study shows that even for relatively weak magnetic fields, the increase in heat transfer for small scale devices is considerable.



---

*Source: 538637-2011-05-23.xml*