Numerical simulation of magnetic nanoparticles targeting in a bifurcation vessel

Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

Contents lists available at ScienceDirect

Journal of Magnetism and Magnetic Materials journal homepage: www.elsevier.com/locate/jmmm

Numerical simulation of magnetic nanoparticles targeting in a bifurcation vessel M.M. Larimi, A. Ramiar n, A.A. Ranjbar Faculty of Mechanical Engineering, Babol Noshirvani University of Technology, PO Box 484, Babol, Iran

art ic l e i nf o

a b s t r a c t

Article history: Received 16 July 2013 Received in revised form 5 October 2013 Available online 18 March 2014

Guiding magnetic iron oxide nanoparticles with the help of an external magnetic field to its target is the principle behind the development of super paramagnetic iron oxide nanoparticles (SPIONs) as novel drug delivery vehicles. The present paper is devoted to study on MDT (Magnetic Drug Targeting) technique by particle tracking in the presence of magnetic field in a bifurcation vessel. The blood flow in bifurcation is considered incompressible, unsteady and Newtonian. The flow analysis applies the time dependent, two dimensional, incompressible Navier–Stokes equations for Newtonian fluids. The Lagrangian particle tracking is performed to estimate particle behavior under influence of imposed magnetic field gradients along the bifurcation. According to the results, the magnetic field increased the volume fraction of particle in target region, but in vessels with high Reynolds number, the efficiency of MDT technique is very low. Also the results showed that in the bifurcation vessels with lower angles, wall shear stress is higher and consequently the risk of the vessel wall rupture increases. & 2014 Elsevier B.V. All rights reserved.

Keywords: Bifurcation Particle tracking Lagrangian method Magnetic nanoparticles Wall shear stress

1. Introduction Fluid flow in bifurcations is a very complex problem of hydromechanics. Unsteady flow and mass transport in bifurcating pipes or vessels are of interest in general and in physiological flows in particular. Flow in bifurcating vessels may be found in blood and lymph-vessels, air-tract and biliary tract. It covers many technical and biomedical applications and many researchers have studied the flow in bifurcation of vessels [1,2] or pipes. Evegren et al. [3], studied Pulsating flow and mass transfer in an asymmetric system of bifurcations. They considered both unsteady- and steady-flow through a three generation system of (non-symmetric) bifurcations. The geometry was consisted of a 901 bifurcation followed by two sets of consecutive symmetric bifurcations. The different inlet conditions affected the flow to the next generation of branches during parts of the cycle. At peak flow and further downstream in the system the effects were negligible. They also found that over a cycle, the mass flow distribution through the outlets can be affected by the inlet velocity conditions and the distribution of a passive scalar is not uniform but depends on the inlet conditions and the Schmidt number. The flow in a carotid bifurcation model was studied by Bharadvaj et al. [4] and Palmen et al. [5]. In most of investigations on flow in large arteries, blood was modeled as a Newtonian fluid.

n

Corresponding author. Tel./fax: þ 98 111 3212268. E-mail address: [email protected] (A. Ramiar).

http://dx.doi.org/10.1016/j.jmmm.2014.03.002 0304-8853/& 2014 Elsevier B.V. All rights reserved.

Gijsen et al. [6] investigated the influence of the non-Newtonian properties of blood on the steady flow in bifurcation vessels. The results showed that in the common carotid artery the nonNewtonian fluid has a flattened axial velocity profile due to its shear thinning behavior. Lee et al. [7] investigated the direct numerical simulation of transitional flow in a carotid bifurcation. They predicted the complex flow field, the turbulence levels and the distribution of biomechanical stresses present in vivo within a carotid bifurcation. Pin et al. [8] simulated the blood flow at arterial bifurcations by the lattice Boltzmann method. Distribution of physical quantities such as the velocity, shear stress and pressure, as well as the location of fluid separation, was investigated numerically. On the other hand the accurate prediction of local dynamical behavior of discrete particles released in the fluid flow is an important key for better understanding and optimization of many processes in numerous branches of science and technology. Examples include a wide spectra of phenomena occurring in environmental (pollution dispersion in atmosphere and oceans), engineering (pharmaceutical industry) and biomedical (deposition of hazardous particles in the human respiratory or cardiovascular system) applications. Zhao et al. [9] adopted discrete trajectory model to simulate particle tracks while the Eulerian method for solving the continuous fluid flow. The results showed that particle deposition and concentration are mainly influenced by the ventilation conditions. It was also found that particles with different sizes (1, 2.5, 5 and10 mm) had different movements in the two ventilated rooms.

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

Nomenclature B P CD Re Cc Rep F T FD u,v FM ! up FB V FG

magnetic magnitude pressure drag coefficient Reynolds number Cunningham correction factor particle Reynolds number total force acting on the particle temperature Drag force velocity components in x and y direction magnetic force velocity of particles Brownian force volume of nanoparticles buoyancy force

0m Fm

χ

Fb

ψ

H

ζ J

ζi KB

υ

M

ρ

mp MnF

59

magnetic permeability of vacuum magnus force magnetic susceptibility Basset force friction coefficient magnetic intensity restitution coefficient current density Gaussian random number Boltzmann constant kinematic viscosity material magnetization density particles mass magnetic numbers

Greek symbols FL

lift force

Longest and Xi [10] studied the effectiveness of direct Lagrangian tracking models for simulating nanoparticle deposition in the upper airways. The objective of their study was to evaluate the effectiveness of direct Lagrangian tracking methods for calculating ultrafine aerosol transport and deposition in flow fields consistent with the upper respiratory tract. Lagrangian deposition results have been compared with a chemical species Eulerian model, which neglects particle inertia. For the tubular and 901 bend geometries, Lagrangian model results with a user defined BM routine agreed well with the Eulerian model. Some new applications stemming from this fusion include pumping and mixing of fluids, as well as the incorporation of switches and valves into lab-on-a-chip devices (Pamme [11]). Magnetic forces are used for transport, positioning, separation and sorting of magnetic as well as non-magnetic objects. Bio-assays have been performed on the surface of magnetic particles trapped inside a micro channel. Many areas in microfluidic applications involve manipulation of particles in a controllable manner. In the design of micro scale devices for cell sorting, cell analysis or cell removal from a sample, it is important to predict the motion of the particles in response to the local flow conditions. For this purpose, various principles and methods have been developed in micro systems, such as the optical tweezers (Furst [12]) and electro kinetic methods (Binyamin et al. [13]). Haverkort et al. [14] performed a Computational simulation of blood flow and magnetic particle motion in a left coronary artery and a carotid artery, using the properties of presently available magnetic carriers. The simulations demonstrated that approximately a quarter of the inserted 4 μm particles can be captured from the bloodstream of the left coronary artery, when the magnet is placed at a distance of 4.25 cm. When the same magnet is placed at a distance of 1 cm from a carotid artery, almost all of the inserted 4 μm particles are captured. There are many researches related to MDT technique for drug delivery in the literatures. Localized medical drug delivery enables a significant local increase of the medical drug in regions affected by disease and leads to a significant reduction of the always present negative side effects of aggressive medical treatments. Experimental studies on animals and preclinical studies on human patients demonstrated potentials of this approach as shown by Lubbe et al. [15] and [16] and Morsi et al. [17]. Riegleret al. [18] presented the

first demonstration of cell targeting using a Magnetic Resonance Imaging (MRI) scanner. They found live human cells, labeled with different iron oxide particles, can be targeted within a vascular vessel model using the magnetic field gradients of an MRI scanner. Liu et al. [19] formulated a hydrodynamic model of Ferro fluids to describe drug carrier flowing in a blood vessel. They found that the enhancement of the magnetic field intensity could slow down the velocity of Ferro fluids and increased the retention of Ferro fluids at the target position. Nacev et al. [20] studied the behavior of ferromagnetic nanoparticles in and around blood vessels under applied magnetic fields. They carried out a detailed analysis to better understand and quantify the behavior of magnetizable particles in-vivo. They found that there are three types of behaviors (velocity dominated, magnetic dominated, and boundary-layer formation) uniquely identified by three essential non-dimensional numbers (the magnetic-Richardson, mass Peclet and Renkin numbers). Kenjeres and Righolt [21] simulated the magnetic capturing of drug carriers in the brain vascular system. They demonstrated that magnetically targeted drug delivery significantly increased the particle capturing efficiency in the target regions. There are two approaches for the numerical simulation of the multiphase flow: Euler–Lagrange and Euler–Euler approach. The Lagrange discrete phase model is based on the Euler–Lagrange approach where the fluid phase is treated as a continuum by solving the time-averaged Navier–Stokes equations, whereas the dispersed phase is solved by numerically integrating the equations of motion for the dispersed phase, i.e. computing the trajectories of a large number of particles or droplets through the calculated flow field. The dispersed phase consists of spherical particles that can exchange mass, momentum and energy with the fluid phase. In general, interactions between particles can be divided in three classes, based on particle volume fraction. According to Elghobashi [22] and Crowe et al. [23], for particle volume fractions less than 10  6, particle motion is influenced by continuous phase properties while practically there is no feedback from the dispersed phase. This class is called “one-way coupling”. For particle volume fractions in the range of 10  6 to approximately 10  3, feedback of the dispersed phase on the properties of the continuous phase fluid dynamics must also be taken into account, this class in known as “two-way coupling”. As the third class, a dense flow is characterized by particle volume fractions

60

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

higher than 10  3. In this condition, which is known as “four-way coupling”, one should also consider particle–particle interactions. In this paper the magnetic particle behavior in an unsteady flow considering the effect of magnetic field in a bifurcation vessel is studied. One way coupling method is used for tracking particles and investigating particle volume fraction in target position and the effect of external magnetic field position. Also the effect of magnetic field on fluid flow and wall shear stresses in bifurcation is studied. All the simulations are performed by the open source CFD software, OpenFoam. The blood flow in the vessel is considered two dimensional, Newtonian, incompressible, laminar and viscose. The forces that effect on the particles are Brownian, Drag, buoyancy force and magnetic force. Partially magnetic field is used in bifurcation walls and its effect is compared for different volume fractions of particles.

magnetic field intensity, H, is given by: qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi γ H ¼ Hx2 þ Hy2 ¼ 2π ðx  aÞ2  ðy bÞ2

The main parameter appearing in the magnetic problems, the magnetic numbers (MnF), due to FHD is considered as: [24]

where ðD=DtÞ is the material derivative, U is the velocity (u, v, 0), P is the pressure, ρ is density and m0 is magnetic permeability of vacuum. The behavior of a bio-magnetic fluid when it is exposed to magnetic field (magnetized) is described by the magnetization property M. Magnetization is the measure of how much the magnetic field affects the magnetic fluid. For isothermal cases, the component of magnetization number is estimated by a linear equation [24]: M ¼ χH

ð2Þ

where χ is a constant, called magnetic susceptibility. Thus, the term ! ðμ0 =ρÞM∇ H in Eq. (1) represents the magnetic force per unit volume and depends on the existence of the magnetic gradient. The governing equations for the x and y components of momentum together with the equation of the continuity in the Cartesian coordinate system under influence of Ferro hydrodynamic (FHD) interactions, may be written as below [24]: ∂u ∂v þ ¼0 ∂x ∂y

ð3Þ

ρ

!
 ∂u ∂u ∂u ∂H þu þv ¼  ∇P þ μð∇2 uÞ þ μ0 M ∂t ∂x ∂y ∂x

ð4Þ

ρ

!
 ∂v ∂v ∂v ∂H þu þv ¼  ∇P þ μð∇2 vÞ þ μ0 M ∂t ∂x ∂y ∂y

ð5Þ

No sleep boundary condition is used for all walls and for the inlet and outlets are used of velocity inlet and pressure outlet boundary condition, respectively. ! ! The terms μ0 Mð∂ H =∂xÞ and μ0 Mð∂ H =∂yÞ in (4) and (5), represent the magnetic force per unit mass along the x and y axes, respectively. The components of the magnetic field intensity Hx and Hy along the x and y coordinates are given, respectively by [24]: Hx ¼

γ xa γ y b ; Hy ¼  2π ðx  aÞ2  ðy bÞ2 2π ðx  aÞ2  ðy  bÞ2

ð6Þ

here (a,b) is the coordinate of the magnetic source position and γ is the magnetic field strength at the point (a,b). The magnitude of the

ð8Þ ð9Þ

where h is the height of the channel. Consequently, for a specific magnetic field, the magnetic numbers can significantly be affected by changing the geometric characteristics and this dependence could be of significant importance from biological view point. According to Eqs. (2) and (8), Mnf can be written as [24]: h μ0 χ H 0 2 h B0 2 χ ¼ 2 v2 ρ v ρμ0 2

ð1Þ

2

B0 ¼ μ0 H 0

Mnf ¼

For two dimensional Newtonian, incompressible, viscous and electrically conducting fluid flow of density ρ (blood) in the presence of uniform transverse magnetic field, created by current carrying wire, the governing equations are as follows: The momentum equations [24]: ! ! μ ! DU 1 ¼  ∇P þ v∇2 U þ 0 M∇ H Dt ρ ρ

h μ0 kH 0 2 h B0 M 0 ¼ v2 ρ v2 ρ 2

MnF ¼

2. Governing equations and boundary conditions 2.1. Fluid phase

ð7Þ

2

ð10Þ

where μ0 is the magnetic permeability of vacuum [4π  10  7 N=A2 ], B0 is the maximum magnetic value and χ is the magnetic susceptibility of fluid in two different conditions; the deoxygenated and oxygenated blood. 2.2. Particles phase In the Lagrangian frame of reference, the equation of motion of a nanoparticle is given by: ! d! up ð11Þ F ¼ dt ! ! where u p is the particle velocity, and F is the acceleration term ! with a unit of force per unit of mass, e.g. N/kg. The F in Eq. (11) is given in the following form which includes contributions from the drag force, magnetic, buoyancy force, Brownian force, thermophoresis force, Saffman's lift force, magnus force and Basset force: F ¼ F D þ F M þ F G þ F B þ F T þF L þF m þ F b

ð12Þ

Most of the above forces should be considered for larger particles and can be neglect for particles in nanoscale. Thermophoretic force is caused by the unequal momentum exchange between the particle and the fluid due to temperature gradient. Since usually human body temperature is about 310 K and the temperature gradient in blood flow is very small, the geometry in this study is considered isothermal and therefore thermophoretic force is neglected. If a particle leads the fluid motion, then the lift force is negative and the particle moves down the velocity gradient towards the wall. Conversely, if the particle lags the fluid, then the lift force is positive and it moves up the velocity gradient away from the wall. This force usually enhances deposition velocity and Wang and Shirazi [25] have shown that neglecting this force has a very small effect in decreasing deposition rate. Intercollision force is the force exerted due to inter-particle collisions. It is usually important when the volume fraction of the dispersed phase is high (more than 10  3). Magnus force is the lift developed due to rotation of the particle. The lift is caused by the pressure difference between both sides of the particle resulting from the velocity difference due to rotation and so, because of neglecting rotation of magnetic particles in blood flow, this force is ignored in this study. The effect of the Basset forces is negligible for particles with density substantially larger than the fluid density. Sommerfeld [26] and Elgobashi and Truesdell [27] showed that the Basset forces are only important for particles with (ρp/ρf⪡1). Thus in this study for particles phase, Magnetic, Drag, buoyancy force and Brownian motion forces are considered as acting forces on particles.

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

2.2.1. Magnetic forces Electromagnetic fields are classically described by Maxwell's equations. For the special case of magneto-static equations that are appropriate for stationary, or slowly varying magnetic fields [18]: ! ! J ¼∇ H ð13Þ ! ∇ Β ¼0

ð14Þ

! ! ! ! ! B ¼ μ0 ð H þ M Þ ¼ μ0 ð H þ χ H Þ

ð15Þ

where B is the magnetic field [T], H is the magnetic intensity [A/m], J is the current density [A/m2], M is the material magnetization [A/m], χ is a dimensionless parameter describing the magnetic susceptibility and μ0 is the permeability of vacuum. A single ferromagnetic particle in a magnetic field will experience a force that depends on the magnetic field and field gradient around it. ! ! ! F m ¼ V  μ0  χ  ð H  ∇ H Þ ð16Þ From Eqs. (11) and (12):

2 ! ∂B ! ! 43  π  a3  χ 6 Bx ∂x ! V χ  ð B ∇Þ B ¼ 6 F m¼ ! 4 μ0 μ0 B By ∂ ∂x

!3 B Bx ∂∂y 7 7 !5 ∂B By ∂y

ð17Þ

where V, is the volume and a is the radius of magnetic particle. 2.2.2. Drag force When the magnetic force of Eq. (17) is applied to a particle, it will make the particle to accelerate in the direction of the ! force until it reaches to an equilibrium velocity, U R , relative to the surrounding fluid. The opposing Stokes drag force on a spherical particle is given by [30]: ! ! ! !  ! 3 ρ f mp ð18Þ C D ð U f  U p Þ U f  U p  F D¼ 4 ρp dp where ρf , ρp , mp and dp are, density of the fluid and particles, ! mass and mean diameter of particles, respectively. U f is the fluid ! velocity and U p is the particle velocity in each time step.C D is the Drag coefficient which can be obtained from the following equation [30]:  24  1 þ 0:15Rep 0:687 CD ¼ ð19Þ Rep where the Reynolds number of particles is defined as [29]: ! !  ρf dp  U f  U p  Rep ¼

μf

ð20Þ

here mf is the fluid viscosity. With combining Eqs. (18)–(20), the drag force acting on a particle will be equal to:  ! ! ! !  ! 3 ρf mp 24    F D¼ ð21Þ  1 þ0:15Rep 0:687 ð U f  U p Þ U f  U p  4 ρp dp Rep 2.2.3. Buoyancy force Buoyancy force is a common issue because the object and surrounding fluid cannot occupy the same physical space simultaneously. Following equation is used to calculate the force on a particle corresponding to buoyancy [30]: ! ρf ! ! F G ¼ mp 1  g ð22Þ

ρp

g is the gravity acceleration, and ρf and ρp are the fluid and particle density, respectively.

61

2.2.4. Brownian motion The components of the Brownian force are modeled as a Gaussian white noise process with spectral intensity Sn;ij [30]: Sn;ij ¼ S0 δi;j where δij is the Kronecker delta function, and ! 216νK B T S0 ¼ π 2 ρdp 2 ðρp =ρÞ2 C c

ð23Þ

ð24Þ

T is the absolute temperature of the fluid, v is the kinematic viscosity, Cc is known as the Cunningham correction factor, KBT is known as the thermal energy (at room temperature¼0.4  10  20 J) and K B is the Boltzmann constant. Amplitudes of the Brownian force components are of the form [29]: rffiffiffiffiffiffiffiffi π S0 ð25Þ F bi ¼ ζ i Δt here ζi is the Gaussian random number with zero mean and unit variance. With combining Eqs. (23)–(25), the Brownian force acting on a particle will be equal to: vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi

u u νK B T uπ 2 216 2 2 t π ρdp ðρp =ρÞ C c F bi ¼ ζ i ð26Þ Δt

2.2.5. Particles motion The particle tracking is performed using a balance between all active forces: ! ! ! ! ! ∑ F ¼ F Dþ F Gþ F Mþ F B

ð27Þ

By considering Eqs. (17), (21), (22) and (25), Eq. (27) can be written as: ! ! ρf ! ! dup ! ! ¼ mp Dð u f  u p Þ þ mp 1  g ∑ F ¼ mp dt ρp rffiffiffiffiffiffiffiffi ! mp  χ ! π S0 ð28Þ ð B  ∇Þ B þ ζ i þ ρp μ0 Δt where: D¼



 3 ρf mp 24v  1 þ 0:15Rep 0:687 4 ρp dp dp

ð29Þ

2.2.6. Collision of particles with the wall Assuming unit vectors n and t to be respectively the normal and tangential unit vectors to the wall, the velocity of the particles after collision with the wall can be determined by: ! ! ! U p ¼ U pn þ U pt

ð30Þ

The normal component of the particle velocity after collision with the wall is evaluated as: ! ! ð31Þ U pn ¼  ς u pn ! u p n is the particle velocity before collision and ς is the wall restitution coefficient where ς A ½0 1. The tangential component of the velocity after the collision with the wall is evaluated as: ! ! ð32Þ U p t ¼ ð1  ψ Þ u p t ! u p t is the tangential velocity before collision and ψ A ½0 1 is the friction coefficient of the wall. The coefficients of restitution and friction are depended on the materials, the surface and the impact velocity and should be determined experimentally. Collision between the particles is neglected in this study.

62

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

3. Geometrical configurations, boundary conditions and model description The geometry used in this study is a bifurcation vessel with one inlet and two outlets, shown in Fig. 1. The height of the inlet and outlets are 1 mm and the length of each branch is considered 1 cm. a uniform velocity was implied for inlet section and a fixed pressure for outlets (pressure outlet B.C). For all walls in both geometry no sleep boundary condition are assumed. Also all walls in the geometry are considered isothermal at 310 K. The magnetic susceptibility of blood (ρ ¼ 1050 ðkg=m3 Þ; ν ¼ 3:1  10  6 m2 s  1 ) is considered 3:5  10  6 for the deoxygenated and  6:6  10  7 for the oxygenated blood, respectively [28,29]. Thus the constant value of χ in the mentioned equations, depended on blood conditions. Although blood is suspension of particles and should be treated as a non-Newtonian flow, it is generally accepted that it behaves as a Newtonian flow in arteries with large diameters (d Z 1 mm) (Huo et al. [31], Nichols and O'Rourke [32]). Furthermore, Pedley [33] and Stroud et al. [34] have shown that in the regions with shear rates higher than 100 s  1 , we can expect Newtonian behavior for the blood flow. Hence in this paper the non-Newtonian behavior of the blood becomes insignificant. In this study, iron oxide particles (Fe3O4) are considered as magnetic particles. Different classes of particles are used, ranging

from 5 nm to 500 nm. The iron oxide particles density (Fe3O4) is assumed to be constant and equal to 6450(kg/m3). A current carrying wire perpendicular to the direction of the flowing fluid is supposed at the position shown in Fig. 1. Horizontal axis (x) in Fig. 1 is the horizontal distance from origin (O).

4. CFD model, mesh and time independency OpenFOAM uses finite volume in a collocated (non-staggered) mesh arrangement, so for solving equation by OpenFoam, the pressure and velocity fields have linked with an algorithm based on Rhie and Chow [35] interpolation. For the stability of the solution, for time scheme, Euler model is adopted, and indicates a first order bounded implicit scheme. For discretization schemes for gradient terms, Gauss linear is applied to all the gradient terms. For the Laplacian operator, the linear interpolation scheme is typically selected. For Time independency, three time steps (5  10  5, 10  6 and 5  10  6) are used and volume fraction in a specified region is compared. The result showed that critically there is no difference in volume fraction, between time step¼10  6 and 5  10  6, and maximum courant number is lower than 1, so for getting this accuracy and satisfying stability criterion, time step¼ 5  10  6 is used. For residual independency checking, four residuals are considered (10  4, 5  10  4, 5  10  5 and 10  6) and maximum shear rate in inner wall

Fig. 1. Geometry of the problem.

Fig. 2. Wall shear stresses in upper wall for Mesh independency investigation at Re¼ 400.

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

are compared and because there were no difference between 5  10  5 and 10  6, the 5  10  5 accuracy is considered. Grid densities of 15,000, 20,000, 27,000 and 33,000 cells were selected to perform a grid independency test as shown in Fig. 2. As a result of occurring recirculation in the upper wall of bifurcation, comparing the Wall Shear Stress (WSS) in this region is selected as a criterion for mesh independency. Comparing the results showed less than 1% difference in WSS between 27,000 and 33,000 grid and hence a 27,000 grid is utilized to the rest of the paper.

63

5. Results and discussion First the effect of the Reynolds number on streamlines and WSS are studied. Fig. 3 shows the streamlines for a bifurcation vessel with an angle of 60□. Increasing the Reynolds number makes the fluid separation to occur at the outer wall of the branch where the secondary flow also occurs. Raising the Reynolds number up to 800 causes the separation zone at the outer walls to become larger and occur earlier.

Fig. 3. Streamlines for various Reynolds numbers at bifurcation with an angle of 60.

Fig. 4. Wall shear stress for upper wall of bifurcation in various Reynolds number.

64

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

Fig. 5. Wall shear stress for inner wall of bifurcation in various Reynolds number.

Fig. 6. Streamlines in various magnetic numbers in deoxygenated (a, b) and oxygenated (c, d) blood. (a) Mnf ¼1.5  105, (b) Mnf ¼ 3.5  105, (c) Mnf ¼  1.5  105 and (d) Mnf ¼  3.5  105.

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

Along the bifurcation, by growing the boundary layer, velocity gradient declines and accordingly WSS decrease. As it can be seen in Fig. 4, there are two points on the upper wall surface where WSS has zero magnitude due to existence of reattachment points. Dimension x in all figures is horizontal distance from origin (O) that is shown in Fig. 1. Fig. 5 and Fig. 6 show the WSS on the upper wall and inner wall along the vessel, respectively. Velocity gradient will increase with increasing the Reynolds number and consequently WSS will increase in the inner and upper wall. In the bifurcation vessels, a stagnation point exists in the middle of the inner wall. Near this point, the pressure and pressure gradient reach to a maximum value along the bifurcation and velocity gradient drop to values near to zero. So, minimum WSS occurs at this point as it is shown in Fig. 5. On the other hand, Fluid flow rate will increase as a reason of forming a recirculation region in separation zone, thus velocity gradient and maximum WSS will happen exactly after the stagnation point on the inner surface (see Fig. 5). Fig. 6 shows the effect of varying magnetic field on the stream lines in oxygenated and deoxygenated blood. By increasing the magnetic number, the number of recirculation areas, backflows and also the length of vorticities will increase. It is also obtained that the effect of the magnetic field on deoxygenated blood is ratherly stronger than the effect of it on oxygenated blood. The ensemble-averaged velocity profiles at the mid-planes of the bifurcation are shown in Fig. 7 at equal magnetic numbers for

65

oxygenated and deoxygenated blood. For both of them, alteration in the direction of velocity vectors in the recirculation zones is depicted in the figure. In the third cross section from left, it is shown that the velocity vectors in the oxygenated blood are in the opposite of the deoxygenated. Also it is revealed that the recirculation zones for the deoxygenated blood take place in the upper branch but for the oxygenated blood, it will happen in the lower branch. The reason lies on the different signs of the magnetic force in the momentum equation for oxygenated and deoxygenated blood. The effect of increasing magnetic number on wall shear stresses at upper wall of bifurcation in oxygenated and deoxygenated blood is depicted in Fig. 8. It is observed in both conditions that with increasing the magnetic number, average WSS in the upper wall will rise because of increasing velocity gradient and existence of secondary flows, as it is shown in Fig. 6. However the effect of magnetic field on blood flow and WSS for deoxygenated blood is stronger than oxygenated blood (Fig. 7b). So it is illustrated that for drug delivery purposes, increasing in external magnetic field has a limit, because of changing in blood flow behavior and recirculation zones formation. These zones will be a big obstacle in magnetic drug targeting. So it is important to use magnetic particles with higher magnetic susceptibility, e.g. iron oxide, due to limits in increasing in external magnetic magnitude. In this paper for magnetic particle targeting, three magnetic magnitude of B¼0.8 T, B¼1.5 T and B¼2.3 T are used and in all conditions, magnetic magnitude don't have any considerable effect on blood flow.

Fig. 7. Velocity vectors in various section for a magnetic number in deoxygenated and oxygenated blood. (a) Deoxygenated and (b) oxygenated.

66

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

Fig. 8. Wall Shear Stresses at the upper wall for deoxygenated blood (a) and average Wall Shear Stress for oxygenated and deoxygenated blood in various magnetic numbers (Mnf) in upper wall.

Fig. 9 shows the effect of various bifurcation angles on inner wall WSS magnitude in the bifurcation. Three angle of bifurcation 301, 601, 901 are considered. As the plot illustrates, decreasing the bifurcation angle intensifies the velocity gradient and consequently, maximum shear stresses. As an example, it is obtained that in comparison with a 301 bifurcation, the maximum WSS in bifurcations with angles of 601 and 901 are decreased 24.86% and 43.75% respectively and maximum shear stress is occurred in bifurcation with angle of 301. So in medical applications, vessels with lower bifurcation angles are more in risk of rupture than vessels with higher angles. But on the other hand by increasing bifurcation angle, maximum pressure near the stagnation point will increase, so the risk of the another medical Phenomenon (Aneurysem), in the part of a blood vessel (artery) will increase or cardiac chamber will swell and either the blood vessel may damage. 5.1. Particle tracing For the first case, three magnetic fields are considered with intensities in the range of 0.8–2.3 T, to deliver particles with 500 nm diameter to a specific region. The target region is considered to be at upper branch. In order to lead the particles to the target, electrical current wire is located at a 1 mm distance from the bifurcation in

the upper branch. The position of electrical current wire and magnetic field lines are shown in Fig. 1. Fig. 10 shows the effect of magnetic field on particles behavior in the blood flow. By influence of magnetic field, most of the nanoparticles are focused on the upper branch. By increasing magnetic field as shown in Fig. 10b and c, number of nanoparticles leaded to the upper branch for delivering to the target region will increase. In the first case, it can be seen that near the recirculation zone at upper branch, nanoparticles follow the blood flow behavior and any of them could not enter to the recirculation zone, but with increasing the magnetic magnitude, particles behavior is changed and nanoparticles volume fraction in the recirculation zone will increase. Red position in bifurcation vessel in Fig. 10 shows the magnitude of magnetic field in this zone. It is obtained that deposition of nanoparticles in the upperwall, near the wire position will increase. So this advantage is used in the next geometry for increasing deposition of nanoparticles in the target region. The variation of volume fraction in the upper and lower branch is exhibited in Table 1. It is observed that by increasing of magnetic magnitude, average volume fraction of particles in upper branch will increase and more percent of injected particles will lead to the upper branch. It is also seen that with increasing magnetic field

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

Fig. 9. Wall shear stress on inner wall for various bifurcation angles.

Fig. 10. Particle behavior under influence of magnetic field in (a) B¼ 0.8 T and (b) B¼ 1.5 T and (c) B ¼2 T.

67

68

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

magnitude from 0.8 T up to 2.3 T, the concentration of particles in the upper branch had a noticeable increasing, and about 90% of injected nanoparticles will be leaded to upper branch. As the second casein this study, geometry shown in Fig. 11, is considered. Target region and wire position for delivering particles to this region are shown in the figure. External magnetic field magnitude and diameter of particles are considered 1.8 T and 500 nm, respectively. All the length and diameter are the same as last geometry. Fig. 12 demonstrates the effect of Reynolds number on magnetic particle targeting. For this purpose, particles in every time steps are injected in blood flow in different Reynolds numbers and the results are compared after 0.5 s. The colored area in Fig. 13 shows magnetic field contours around the current carrying wire. As shown in Fig. 12, increasing the Reynolds number causes the volume fraction of magnetic particles in target region to reduce. Table 2 shows the average volume fraction in the target region for various Reynolds numbers. Obviously, average volume fraction

and efficiency of magnetic particle targeting will decrease in target region at higher Reynolds numbers. As it is shown in Table 2, raising the Reynolds number from Re ¼150 up to Re¼550, causes the average volume fraction to decrease more than 800%. So, it is illustrated that MDT technique is not an effective medicine method in vessels with higher Reynolds numbers. maximum pick for Re¼150 and Re¼400, in Fig. 12, shows that the deposition of nanoparticles in this zone, which are more affected by the magnetic field, is higher than other areas and in the medical applications these zones can be areas that are affected by cancers or tumors. Fig. 13 shows the effect of particle's diameter on delivering magnetic particles to the target region. For the first case (a), distribution of magnetic particles with 500 nm diameter at Re¼150 and B¼1.8 T in deoxygenated blood, is investigated. It can be observed that some of the injected magnetic particles will be delivered to the target region by influence of external magnetic field. In second case (b), particles with 50 nm diameter are checked out and found that

Table 1 Average nanoparticles volume fraction in various magnetic field in bifurcation. Magnetic magnitude

0.8 T

Volume fraction of particles in upper branch Volume fraction of particles in lower branch Percent of concentration of nanoparticles in the upper branch to that of lower branch

1.5 T 6

9.859361  10 8.170779  10  6 20.5%

Fig. 11. Second geometry.

Fig. 12. Volume fraction in target region for various Reynolds number.

2.3 T 5

1.360218  10 6.163491  10  6 105%

1.783298  10  5 2.158773  10  6 720%

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

69

Fig. 13. Magnetic particle distribution near target region for particles with various diameter.

Table 2 Average volume fraction in target region at various Reynolds numbers. Reynolds number

Re¼150

Re¼ 400

Re¼ 550

Average volume fraction in target region

1.942496  10  5

7.439026  10  6

1.955779  10  6

number of particles delivered to target region is decreased. In the last case (c), the diameter of particles is reduced to 5 nm. It is observed that at this Reynolds and magnetization number, none of injected particles are delivered to the target region. The reason is that, by reducing the diameter of magnetic particles, in Eq. (13), volume of nanoparticles is decreased and consequently magnetic force that effect on each particle is reduced and nanoparticles cannot be lead to the secondary branch. In the last part, the effect of ratio of particle diameter (Dp) to magnetic core diameter (Dm) on delivery efficiency is demonstrated (Fig. 14). In drug delivery applications as it was mentioned, the magnetic cores are coated with a layer of drug and certain biocompatible polymers such as dextran or polyethylene glycol. Therefore, it is always necessary to have Dp/Dm 41 in order to be able to deliver some amount of the medical drug attached to the magnetic carrier. In Figs. 14 and 15, for Dp/Dm ¼1, 2,3,4,5, 10, and 1, particles volume fraction in target region after 0.5 s is obtained. For the first case (Dp/Dm ¼ 1, Fig. 14a), none of the injected particles, are delivered to the target region. With increasing the diameter

ratio, the efficiency of drug delivery will increase. The effect of Dp/Dm on volume fraction is illustrated in Fig. 15. As it is observed, raising the diameter ratio reduces the volume fraction. Decreasing the magnetic core, declines the magnetic force magnitude that affect on particles, because the volume of magnetic core in Eq. (15) will decrease but the drag force that act on the outer diameter of particles (coated layer) will be constant. So the magnetic force magnitude is not enough to construct to drag force and consequently MDT technique efficiency will decrease.

6. Conclusions In this study the behavior of blood flow and magnetic particles (Fe3O4) influenced by external magnetic field in a bifurcation vessel is investigated. The flow is incompressible, unsteady, Newtonian, viscous and electrical conductive. It is observed that in high Reynolds numbers, secondary flow and separation occurs and by increasing the Reynolds number as it showed in Fig. 3,

70

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

Fig. 14. Magnetic particle distribution near target region for various Dp/Dm.

Fig. 15. Particle average volume fraction versus (Dp/Dm).

length of the recirculation zone increases and back flows occurs earlier. Also effect of high magnetization number on blood flow in two cases, oxygenated blood and deoxygenated blood is investigated. It is observed that by increasing external magnetic field (increasing electrical current in wire) recirculation and back flows and their length are increased in both cases, but deoxygenated blood, is more affected by magnetization number in comparison with oxygenated blood. For the particle phase, effect of magnetic field gradient on behavior of magnetic nanoparticles

is investigated by use of a Lagrangian particle tracking model. It is illustrated that, by use of an external magnetic field, particles can be delivered to a target region, and with increasing magnetic gradient, volume fraction of particles that are delivered to the target region, is increased. By increasing the Reynolds number, volume fraction of nanoparticles delivered to targeted region will be reduced. Also it is shown that by increasing magnetic core in drug delivery technique, the efficiency of MDT technique increases.

M.M. Larimi et al. / Journal of Magnetism and Magnetic Materials 362 (2014) 58–71

References [1] F. Ghalichi, X. Deng, A.D. Champlain, Y. Douville, M. King, R. Guidoin, Low Reynolds number turbulence modeling of blood flow in arterial stenosis, Biorheology 35 (1998) 281–294. [2] J.S. Stroud, S.A. Berger, D. Saloner, Numerical analysis of flow through a severely stenotic carotid artery bifurcation, J. Biomech. Eng. 124 (2002) 9–20. [3] P. Evegren, J. Revstedt, L. Fuchs, Pulsating flow and mass transfer in an asymmetric system of bifurcations, Comput. Fluids 49 (2011) 46–61. [4] B. Bharadvaj, R. Mabon, D. Giddens, Steady flow in a model of the human carotid bifurcation: Part I—Flow visualization, J. Biomech. Eng. 15 (1982) 349–362. [5] D.E.M. Palmen, F.J.H. Gijsen, F.N.V. Vosseand, J.D. Janssen, Diagnosing minor stenosis in carotid artery bifurcation models using the disturbed velocity field, J. Vasc. Surg. 3 (1997) 26–41. [6] F.J.H. Gijsen, F.N. van de Vosse, J.D. Janssen, The influence of the nonNewtonian properties of blood on the flow in large arteries: steady flow in a carotid bifurcation model, J. Biomech. 32 (1998) 601–608. [7] S.E. Lee, S.W. Lee, P.F. Fischer, H.S. Bassiounyd, F. Lothe, Direct numerical simulation of transitional flow in a stenosed carotid bifurcation, J. Biomech. 41 (2008) 2551–2561. [8] J.Y. Pin, K.X. Ying, L.D. He., The blood flow at arterial bifurcations simulated by the lattice Boltzmann method, Chin. Phys. Lett. 7 (2009) 074702. [9] B. Zhao, Y. Zhang, X. Li, X.g Yang, D. Huang, Comparison of indoor aerosol particle concentration and deposition in different ventilated rooms by numerical method, Build. Environ. 39 (2004) 1–8. [10] P.W. Longest, J. Xi, Effectiveness of direct Lagrangian tracking models for simulating nanoparticle deposition in the upper airways, Aerosol Sci. Technol. 41 (2007) 4. [11] N. Pamme, Magnetism and microfluidics, Lab Chip 6 (2006) 24–38. [12] E.M. Furst, Applications of laser tweezers in complex fluid rheology, Adv. Colloid Interface Sci. 10 (2005) 79–86. [13] G. Binyamin, T.D. Boone, H.S. Lackritz, A.J. Ricco, A.P. Sassi, S.J. Williams, R.E. Oosterbroek, B.A. van den, Plastic microfluidic devices, Electro kinetic manipulations, life science applications, and production technologies, Lab Chip (2003) 83–112. [14] D.J.W. Haverkort, S. Kenjeres, C.R. Kleijn, Computational simulations of magnetic particle capture in arterial flows, Ann. Biomed. Eng. 37 (2009) 2436–2448. [15] A.S. Lubbe, C. Bergemann, H. Riess, F. Schriever, P. Reichardt, Clinical experiences with magnetic drug targeting: a phase study with 40 -epidoxorubicin in 14 patients with advanced solid tumors, Cancer Res. 56 (20) (1996) 4686–4693. [16] A.S. Lubbe, C. Bergemann, W. Huhnt, T. Fricke, H. Riess, J.W. Brock, D. Huhn, Preclinical experiences with magnetic drug targeting: tolerance and efficacy, Cancer Res. 56 (20) (1996) 4694–4701.

71

[17] S.A. Morsi, A.J. Alexander, An investigation of particle trajectories in two-phase flow systems, J. Fluid. Mech. 55 (1972) 193–208. [18] J. Riegler, J.A. Wells, P.G. Kyrtatos, A.N. Price, Q. Pankhurst, M. Lythgoe, Targeted magnetic delivery and tracking of cells using a magnetic resonance imaging system, Biomaterials 31 (2010) 5366–5371. [19] H.D. Liu, W. Xu, S.G. Wang, Z.J. Ke, Hydrodynamic modeling of ferro fluid flow in magnetic targeting, Appl. Math. Mech. 29 (2008) 1341–1349. [20] A. Nacev, C. Beni, O. Bruno, B. Shapiro, The behaviors of ferromagnetic nanoparticles in and around blood vessels under applied magnetic fields, J. Magn. Magn. Mater. 323 (2011) 651–668. [21] S. Kenjeres, B.W. Righolt, Simulations of magnetic capturing of drug carriers in the brain vascular system, Int. J. Heat Fluid Flow 35 (2012) 68–75. [22] S. Elghobashi, On predicting particle-laden turbulence flows, Appl. Sci. Res. 52 (1994) 309–329. [23] C. Crowe, T. Troutt, J. Chung, Numerical models for two-phase turbulent flows, Annu. Rev. Fluid Mech. 28 (1996) 11–43. [24] E.E. Tzirtzilakis, V.D. Sakalis, N.G. Kafoussias, P.M. Hatzikonstantinou, Biomagnetic fluid flow in a 3D rectangular duct, Int. J. Numer. Methods Fluids 44 (2004) 1279–1298. [25] J. Wang, S.A. Shirazi, ACFD Based, Correlation for erosion factor for long-radius elbows and bends, J. Energy Res. Technol. 125 (2003) 26–34. [26] M. Sommerfeld, Theoretical and Experimental Modeling of Particulate Flow, VKI Lecture Series (2000) 06. [27] S. Elgobashi, G.C. Truesdell, Direct simulation of particle dispersion in decaying isotropic turbulence, J. Fluid Mech. 242 (1992) 655–700. [28] M. Motta, Y. Haik, A. Gandhari, C.J. Chen, High magnetic field effects on human de-oxygenated hemoglobin light absorption, Bioeletrochem. Bioenerg. 47 (1998) 297. [29] Y. Haik, V. Pai, C.J. Chen, Biomagnetic fluid dynamics, in: W. Shyy, R. Narayanan (Eds.), Fluid Dynamics at Interfaces, Cambridge University Press, 1999, pp. 439–452. [30] D. Wen, L. Zhang., Y. He, Flow and migration of nanoparticle in a single channel, Heat Mass Transfer 45 (2009) 1061–1067. [31] Y. Huo, J.S. Choy, M. Svendsen, A.K. Sinha, G.S. Kassab, Effects of vessel compliance on flow pattern in porcine epicardial right coronary arterial tree, J. Biomech. 42 (2009) 594–602. [32] W.W. Nichols, M.F. O'Rourke, McDonald's Blood Flow in Arteries: Theoretical, Experimental and Clinical Principles, sixth ed., CRC Press, 2011. [33] T.J. Pedley, The Fluid Mechanics of Large Blood Vessels, Cambridge University Press, 1980. [34] J.S. Stroud, S.A. Berger, D. Saloner, Influence of stenosis morphology on flow through severely stenotic vessels: implications for plaque rupture, J. Biomech. 33 (2000) 443–455. [35] C.M. Rhie, W.L. Chow, Numerical study of the turbulent flow past an airfoil with trailing edge separation, AIAA J. 21 (1983) 1525–1532.