Coupled field analysis
The dynamic characteristics of magnetic particles are mainly determined by two factors in microfluidic chips, the gradient magnetic field applied externally and the flow field in the microchannel. Therefore, the motion characteristics of magnetic beads in the microchannel mainly involve the theory of magnetostatics, hydrodynamics and Newton’s second law.
Magnetic field
Magnetic field is the main way to control material transport in microfluidic chip. Compared with electric field, magnetic field has unique advantages: magnetic force is not affected by sample concentration, pH value and other parameters; No heat is generated by permanent magnet; Magnetic field can penetrate most of the materials such as glass and polymer36; magnetic field is mild and non-destructive to most living cells. In magnetic beads preset technology, compound magnetic field effect of different permanent magnets on the magnetic beads is mainly considered.
Flow field
When the magnetic beads enter the chip from outside, the capture efficiency of magnetic beads is closely related to the flow rate. When the flow rate is low, the capture efficiency is high, but the efficiency of biological analysis system is low. When the flow rate is high, the capture efficiency is low. Although the suitable flow rate and capture efficiency can be achieved with optimizing, using the magnetic beads preset technology can weaken the influence of flow rate, focus on the influence of static flow field on magnetic beads, and simplify the analysis model.
Mathematical model of magnetic beads motion
The external forces on magnetic beads in the microchannel include magnetic force, fluid viscous force, interaction force between particles, gravity and buoyancy37,38. According to Newton’s second law:
$$ mfrac{{d^{2} tau }}{{dt^{2} }} = F_{m} + F_{d} + F_{g} + F_{f} $$
(1)
where m is the mass of particle, τ is the position vector diameter of particle, and Fm, Fd, Fg, Ff are the magnetic force, viscous force, gravity and buoyancy. Since magnetic force and viscous force are much greater than gravity and buoyancy, the model is simplified as formula (2), ignoring gravity and buoyancy.
$$ mfrac{{d^{2} tau }}{{dt^{2} }} = F_{m} + F_{d} $$
(2)
The forces on magnetic bead are shown in Fig. 9. The magnetic field force is decomposed into three components (Fx, Fy, Fz).
Magnetic force F
m
Magnetic force Fm can be expressed as:
$$ F_{m} = mu_{0} MV_{p} nabla H $$
(3)
where μ0 is the vacuum permeability, M is the field dependent magnetization of magnetic bead, Vp is the volume of magnetic bead, and H is the space magnetic field strength. When the particle is saturated magnetized, its field dependent magnetization is the saturation magnetization, that is, M = Ms.
In general, the field dependent magnetization follows the formula:
where H is the mode of magnetic field intensity, H = |h|, and f(H) satisfies the equation:
$$ f(H) = left{ {begin{array}{*{20}l} {frac{3chi }{{3 + chi }},} hfill & {H > frac{chi + 3}{{3chi }}M} hfill \ {frac{{M_{s} }}{H} = chi ,} hfill & {H ge frac{chi + 3}{{3chi }}M} hfill \ end{array} } right. $$
(5)
where, (V_{p} = frac{3}{4}pi r^{3}), (chi) is the magnetic susceptibility of the magnetic bead, and R is the radius of the magnetic bead.
Because (H = frac{B}{{mu_{0} }}), so:
$$ F_{m} = frac{1}{{mu_{0} }}f(H)V_{p} Bnabla B{ = }frac{1}{{mu_{0} }}f(H)V_{p} left( {B_{x} frac{partial B}{{partial x}} + B_{y} frac{partial B}{{partial y}} + B_{z} frac{partial B}{{partial z}}} right) $$
(6)
where B is the magnetic flux density.
The magnetic field force on magnetic bead depends not only on the magnetic induction intensity, but also on the magnetic field gradient. Magnetic bead tends to move towards the maximum magnetic field intensity in the non-uniform magnetic field. The magnetic force in different directions are as follows:
$$ left{ begin{gathered} F_{{text{x}}} = frac{1}{{mu_{0} }}f(H)V_{p} Bnabla B{ = }frac{1}{{mu_{0} }}f(H)V_{p} left( {B_{x} frac{{partial B_{x} }}{partial x} + B_{y} frac{{partial B_{x} }}{partial y} + B_{z} frac{{partial B_{x} }}{partial z}} right) hfill \ F_{y} = frac{1}{{mu_{0} }}f(H)V_{p} Bnabla B{ = }frac{1}{{mu_{0} }}f(H)V_{p} left( {B_{x} frac{{partial B_{y} }}{partial x} + B_{y} frac{{partial B_{y} }}{partial y} + B_{z} frac{{partial B_{y} }}{partial z}} right) hfill \ F_{{text{z}}} = frac{1}{{mu_{0} }}f(H)V_{p} Bnabla B{ = }frac{1}{{mu_{0} }}f(H)V_{p} left( {B_{x} frac{{partial B_{z} }}{partial x} + B_{y} frac{{partial B_{z} }}{partial y} + B_{z} frac{{partial B_{z} }}{partial z}} right) hfill \ end{gathered} right. $$
(7)
Viscous force F
d
In microfluidic system with low Reynolds number, the viscous force acting on particles can be obtained by Stokes law:
$$ F_{d} = 6pi eta r(v_{f} – v_{p} )f_{D} $$
(8)
where, η is the dynamic viscosity of fluid, r is the radius of magnetic bead, vf is the velocity of fluid, vp is the velocity of magnetic bead, and fd is the hydrodynamic resistance coefficient.
$$ v_{p} = frac{dtau }{{dt}} $$
(9)
When the wall effect is considered, the expression of hydrodynamic resistance coefficient is as follows:
$$ f_{D} = left[ {1 – frac{9}{16}left( {frac{r}{r + s}} right) + frac{1}{8}left( {frac{r}{r + s}} right)^{3} – frac{45}{{256}}left( {frac{r}{r + s}} right)^{4} – frac{1}{16}left( {frac{r}{r + s}} right)^{5} } right]^{ – 1} $$
(10)
where s is the distance between the magnetic bead and the wall of microchannel. When the magnetic bead is far away from the wall, fD = 1; When the magnetic bead is near the wall, fD is slightly greater than 1.
Speed of magnetic bead
The time required for magnetic bead to reach a new equilibrium is extremely small. It can be considered that the motion of magnetic bead in the microchannel is always in a Quasi equilibrium state. The motion equation of the magnetic bead can be simplified as:
$$ F_{m} + F_{d} { = }0 $$
(11)
By introducing Eqs. (8) and (9), we get the following results:
$$ v_{p} = frac{dtau }{{dt}} = frac{1}{{6pi eta rf_{D} }}F_{m} + v_{f} $$
(12)
The velocity components in each direction can be expressed as:
$$ left{ begin{gathered} v_{px} = frac{dx}{{dt}} = frac{1}{{6pi eta rf_{D} }}F_{mx} + v_{fx} hfill \ v_{py} = frac{dy}{{dt}} = frac{1}{{6pi eta rf_{D} }}F_{my} + v_{fy} hfill \ v_{pz} = frac{dz}{{dt}} = frac{1}{{6pi eta rf_{D} }}F_{mz} + v_{fz} hfill \ end{gathered} right. $$
(13)
When the magnetic beads with a certain velocity flows go through the microchannel, the magnetic beads will be shifted and separated from the fluid by magnetic force. Considering that the width to height ratio of the microchannel is large, the separation process can be simplified as a 2D model.
Motion analysis of preset magnetic bead capture
The capture motion includes magnetic beads captured by permanent magnet 1 (preset permanent magnet) after magnetic beads released and fixed on the upper surface of the channel during solution washing. The force analysis of magnetic bead capture is shown in Fig. 10a, where Ff is the friction force between magnetic bead and the channel upper surface. When the magnetic bead is stable in the capture position, Fmx1 = 0 and vP = 0. The condition for the successful capture of magnetic beads under the velocity vf is as follows:
$$ F_{f} { = }F_{{text{d}}} $$
(14)
Force analysis of preset magnetic beads. (a) Preset magnetic bead capture. (b) Upward mixing motion. (c) Downward mixing motion.
In the critical state:
$$ F_{f} = fF_{{{text{mz1}}}} $$
(15)
where f is the friction coefficient between the magnetic bead and the flow channel. Introduce Eqs. (8) and (15) into Eq. (14), and get:
$$ F_{mz1} { = }frac{{6pi eta rf_{D} }}{f}v_{f} $$
(16)
Let (frac{{6pi eta rf_{D} }}{f} = [{text{a}}]), When the microfluidic system is determined, [a] is a certain value, then:
$$ F_{mz1} = [a]v_{f} $$
(17)
Therefore, only when the z-direction magnetic force produced by permanent magnet 1 is greater than [a] times of the liquid flow rate, the magnetic beads can be effectively captured.
Motion analysis of preset magnetic bead mixing
The mixing motion between magnetic beads and solution can be divided into upward mixing motion and downward mixing motion, which are driven by different permanent magnets. Because of using preset magnetic beads technology, vfz is considered as 0.
Upward mixing motion
The upward mixing motion is achieved by permanent magnet 1. At this time, the permanent magnet 2 falls far enough, so the influence of its magnetic field on magnetic bead can be ignored. The force analysis of upward mixing motion of magnetic bead is shown in Fig. 10b. The achievement of upward mixing motion should satisfy Fmz1 > Fd, considering the quasi equilibrium state of magnetic bead motion:
$$ v_{pz1} = frac{dz}{{dt}} = frac{1}{{6pi eta rf_{D} }}F_{mz1} $$
(18)
The microchannel height is Q, and integrate both sides:
$$ t_{1} = frac{{6pi eta rf_{D} Q}}{{F_{mz1} }} $$
(19)
where t1 is the time required for upward mixing motion.
Downward mixing motion
The downward mixing motion is achieved by permanent magnet 2. The force analysis of downward mixing motion of magnetic bead under the influence of the magnetic field of permanent magnet 1 is shown in Fig. 10c. The achievement of downward mixing motion should satisfy fmz2 > fmz1 + Fd, considering the quasi equilibrium state of magnetic bead motion:
$$ left{ begin{gathered} v_{pz1} = frac{1}{{6pi eta rf_{D} }}F_{mz1} hfill \ v_{pz2} = frac{1}{{6pi eta rf_{D} }}F_{mz2} hfill \ end{gathered} right. $$
(20)
The downward velocity of magnetic bead is:
$$ v_{pz} = frac{dz}{{dt}}{ = }v_{pz2} – v_{pz1} = frac{1}{{6pi eta rf_{D} }}(F_{mz2} – F_{mz1} ) $$
(21)
Integrate both sides:
$$ t_{2} = frac{{6pi eta rf_{D} Q}}{{F_{mz2} { – }F_{mz1} }} $$
(22)
where t2 is the time required for downward mixing motion.
The mixing efficiency of magnetic beads is determined by t1 and t2. After the magnetic beads are constant, the mixing efficiency depends on the following factors: (1) the smaller fluid dynamic viscosity, the higher mixing efficiency; (2) the greater magnetic field intensity of preset permanent magnet, the shorter upward mixing motion time, but the downward mixing motion needs to be considered comprehensively; (3) the greater magnetic field intensity difference, the shorter downward mixing motion time; (4) the smaller channel height, the higher mixing efficiency. After the type of magnetic beads, the size and quantity of flow channel and the type of different reaction solutions are determined, the specification of permanent magnet should be reasonably selected to ensure the magnetic beads motion effect, simplify the device structure and reduce the structure size.
Through the derivation of the upward mixing movement time and downward mixing movement time of preset magnetic beads, the time required for a single full mixing movement can be clearly obtained, and the total time required for the mixing movement can be obtained according to the number of up-and-down movements through related DNA extraction and elution effect experiment, which provides a basis for the determination of relevant parameters of the control magnetic field.
Numerical simulation for magnetic field
In the magnetic beads preset technology, the most important is how to select the appropriate permanent magnet to ensure the capture efficiency and mixing efficiency. We have developed a method to simulate the magnetic field by using COMSOL and verify the variation characteristics of magnetic field by changing the aspect ratio of permanent magnet in the microfluidic channel. As shown in Fig. 11a, the 2D model was composed of permanent magnet, microchannel, chip and film. The related dimensions are: the length of microchannel 2 mm, the width of microchannel 0.1 mm, the permanent magnet initial dimension of width w = 0.15 mm and height h = 0.15 mm. The material of microfluidic channel was set as air, which could truly simulate the impact of the environment. The boundary condition of the air was set as magnetic insulation, and the grid was divided according to the standard. The relevant parameters are: relative permeability of permanent magnet 1, vacuum permeability 4π × 10−7 H A−2, permanent magnet x axial remanence 0, permanent magnet y axial remanence ± 0.4 T.
Numerical simulation for magnetic field. (a) 2D analysis model. (b) Spatial distribution of magnetic field. (c) Magnetic field intensity component in microchannel centerline. (d) Magnetic field x components of microchannel centerline with different aspect ratio.
The spatial distribution of magnetic field is shown in Fig. 11b. The densest magnetic flux density region is concentrated in the interior of permanent magnet and gradually decreases upward along the radial direction of the channel. The gradient magnetic field in the microchannel above the permanent magnet is larger and the effect of the magnetic field extending to both sides is reduced. Therefore, the magnetic flux density near the magnet direction in the microchannel is larger and the magnetic field is stronger.
The magnetic field intensity component generated by permanent magnet in microchannel centerline is shown in Fig. 11c. The magnetic field gradient in x-direction is the largest near x = 1 mm, and the magnetic field intensity in y-direction reaches the maximum. Since the magnetic field force is directly proportional to the magnetic field gradient, it can be inferred that the magnetic force in x-direction will have the maximum near x = 1 mm. Therefore, the magnetic beads will tend to move in the central region of the permanent magnet.
Magnetic field x components of microchannel centerline with different aspect ratio are shown in Fig. 11d. When the width w of the permanent magnet remains constant, the aspect ratio has no effect on the curve shape of the magnetic field intensity, which is a sinusoidal curve. The maximum value of x-direction magnetic field intensity and gradient magnetic field increase with the increase of magnet height h, but when the height of permanent magnet increases to a certain value, the increase of magnetic field intensity in the channel decreases gradually until it keeps constant. Therefore, the aspect ratio of permanent magnet can be increased appropriately to obtain greater magnetic field force and improve the magnetic beads capture and mixing efficiency.
