Next Article in Journal
Direct Shear Stress Mapping Using a Gallium Nitride LED-Based Tactile Sensor
Next Article in Special Issue
Continuous On-Chip Cell Washing Using Viscoelastic Microfluidics
Previous Article in Journal
Stability Compensation Design and Analysis of a Piezoelectric Ceramic Driver with an Emitter Follower Stage
Previous Article in Special Issue
A Method for Rapid, Quantitative Evaluation of Particle Sorting in Microfluidics Using Basic Cytometry Equipment
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Numerical Study of Viscoelastic Microfluidic Particle Manipulation in a Microchannel with Asymmetrical Expansions

1
Department of Hydraulic Engineering, College of Civil Engineering and Architecture, Zhejiang University, Hangzhou 310058, China
2
School of Mechanical and Mining Engineering, The University of Queensland, St Lucia, QLD 4072, Australia
3
School of Electrical and Electronic Engineering, Nanyang Technological University, Singapore 639798, Singapore
*
Authors to whom correspondence should be addressed.
Micromachines 2023, 14(5), 915; https://doi.org/10.3390/mi14050915
Submission received: 29 March 2023 / Revised: 20 April 2023 / Accepted: 21 April 2023 / Published: 23 April 2023
(This article belongs to the Special Issue Viscoelastic Microfluidics and Cell Sorting)

Abstract

:
Microfluidic microparticle manipulation is currently widely used in environmental, bio-chemical, and medical applications. Previously we proposed a straight microchannel with additional triangular cavity arrays to manipulate microparticles with inertial microfluidic forces, and experimentally explored the performances within different viscoelastic fluids. However, the mechanism remained poorly understood, which limited the exploration of the optimal design and standard operation strategies. In this study, we built a simple but robust numerical model to reveal the mechanisms of microparticle lateral migration in such microchannels. The numerical model was validated by our experimental results with good agreement. Furthermore, the force fields under different viscoelastic fluids and flow rates were carried out for quantitative analysis. The mechanism of microparticle lateral migration was revealed and is discussed regarding the dominant microfluidic forces, including drag force, inertial lift force, and elastic force. The findings of this study can help to better understand the different performances of microparticle migration under different fluid environments and complex boundary conditions.

1. Introduction

Microfluidic particle manipulation technologies are widely used in many areas, including biomedical engineering, material sciences, food engineering, and pharmaceuticals [1,2,3,4]. Generally, the microfluidic manipulation methods can be roughly classified into two categories, active methods and passive methods [5,6].
For active methods, different types of force fields can be used for particle manipulation, such as electric fields [7,8], magnetic fields [9,10], and acoustic waves [11,12]. However, the difficulties of accurate integration and maintenance involving microfluidic channels and the additional devices that apply the force fields have been found to be troublesome in different applications [13,14,15]. By comparison, passive methods only rely on the hydrodynamics that are generated in the fluid domain to realize the microparticle manipulation, which requires a deep understanding of the physics of fluids and fluid–particle interactions. Well-designed microchannels of specific geometric characteristics can introduce an expected hydraulic force field to manipulate microparticles without assistance from additional physical fields [16,17,18]. Lee et al. [19] summarized a review of the passive methods of microfluidic mixing, introducing the structural geometrics of a lamination-based micromixer [20], T-shaped micromixer [21], Planar Asymmetric Split-and-Recombine micromixer [22], etc. Furthermore, the tasks of differentiating, arranging, and isolating varying microparticles can present a greater level of complexity, while also offering a wider range of potential uses. For instance, Preira et al. [23] proposed a microfluidic method for the passive sorting and separation of non-adherent cell populations by deformability. The method can be adapted in hospitals with only the microchip and some standard syringes and pumps. Tan et al. [24] proposed a design of a passive microfluidic channel to realize the control of droplet volume and sorting. The droplet sorting was achieved using a bifurcating flow design that enables the separation of droplets according to their size.
Viscoelastic microfluidics has emerged as an exciting and rapidly growing field in recent years because the intrinsic properties of viscoelastic fluids can lead to unique particle migration phenomena. The growing interest in this area includes developing particle manipulation techniques and applications of viscoelastic microfluidics [25,26,27]. For example, Zhang et al. [28] successfully isolated drug-treated Escherichia coli by shape-based separation using viscoelastic microfluidics. Their study confirmed that two main factors, namely, the sheath-to-sample flow rate ratio and the concentration of polyethylene oxide (PEO), had a significant effect on the separation effectiveness of the viscoelastic microfluidic device. Feng et al. [29] investigated particle focusing and separation in a spiral channel and explained the variation in particle focusing position based on the coeffects of inertial flow, viscoelastic flow, and Dean flow.
Less fundamental and simulation research is undertaken in viscoelastic microfluidics, which is essential for the better understanding of the behavior of these fluids and particles, and the prediction of their flow behavior in microscale systems, for both fundamental scientific advances and practical applications. The numerical models used for viscoelastic fluid flow, such as the Oldroyd-B, Giesekus, and Phan-Thien/Tanner models [30,31,32,33], are valuable tools that help explain the complex effects of fluid rheology on particle migration. In terms of applications, most simulations are performed in straight microchannels with square cross-sections [33,34,35] to investigate the intricate effects of elastic and inertial forces. In addition, most of these papers focus on a low Reynolds number and Weissenberg number. Di et al. [36] explored particle elasto-inertial focusing in straight microchannels through numerical simulation. The study scrutinized particle focusing under multiple control parameters, such as the Reynolds number, Weissenberg numbers, and particle diameter, to shed light on the underlying force competition mechanisms. Mohammad et al. [37] employed the Giesekus constitutive equation to examine the impact of flowrate, particle size, and shear-thinning degree on the resulting focusing patterns. To accomplish this, experiments and 3D simulations were conducted. The results reveal that, at low flowrates, particles were guided towards the channel center by the combined forces of elastic force and secondary flow. Conversely, at high flowrates, the increased shear-gradient lift and reversed elastic force direction caused particles to be dispersed away from the center.
While considerable progress has been made through experimentation in expanded and contracted cavity array (ECCA) channels, these microchannels have not been given the same level of attention in numerical simulations. Consequently, there remains a lack of understanding regarding the mechanisms and interrelationships between different types of microfluidics, and a thorough and systematic examination is therefore warranted. In our previous study [38], we explored the particle migration behavior in viscoelastic fluids in a straight microchannel with additional triangular cavity arrays. The experimental results indicated that the microparticle trajectory is sensitive to both the flow rate and fluid viscoelasticity. However, the mechanism still remained unclear due to the limitation of sight views in experiment.
In this study, we introduce numerical simulation to provide a global transparent view of different scenarios of the fluid domain characteristics and the microparticles’ trajectory in the aforementioned microchannels. The numerical model was established in COMSOL Multiphysics 6.0 and validated by the experimental results. According to the simulation results, the effects of flow rate and viscoelasticity on the introduced hydraulic forces fields were analyzed and are discussed in detail. The root reason for the lateral migration of microparticles and its intensity was revealed by extracting and quantitatively comparing the dominant hydraulic forces, including the drag force, the inertial lift force, and the elastic force. Based on the proposed numerical simulation methods and the findings of this study, one can undertake the optimization of the structural design of the microfluidic channels, identification of the optimal flow rate, and determination of the fluid physical characteristics such as viscoelasticity.

2. Mathematical Model

2.1. Model Design and Establishment

The schematic figure of the microchannel used in this study is shown in Figure 1. The microchannel domain is a straight channel with a cross-section of 100 μm × 40 μm (width × height), and 26 additional cavities on one side with the shape of right-angled isosceles triangles. The longest side of the triangular cavities is 900 μm, the spacing between two adjacent cavities is 900 μm, and the full length of the straight channel is 48 mm. Since the channel height is much less than the length and the width (H << L, H << W), the model can be simplified to a 2D solution. Details of the chip fabrication can be found in the reference [38].
COMSOL Multiphysics 6.0 was used to establish this corresponding numerical simulation model. The stationary study of the fluid field was undertaken using the viscoelastic fluid module. Since the Reynolds number in the fluid domain is much less than 2000, the fluid conditions were considered as stationary laminar flow. Subsequently, the results of the stationary calculations were introduced into the particle-tracing module to calculate the real-time positions and velocities of the migrating suspended microparticles at different time steps.
In previous work [38], we experimentally investigated continuous plasma extraction in the ECCA channel using a blood sample diluted with 500 ppm, 1000 ppm, and 2000 ppm polyethylene oxide (PEO) solutions. The estimated viscosity values for 500 ppm, 1000 ppm, and 2000 ppm PEO solutions were 2 × 10−3, 3 × 10−3, and 5 × 10−3 Pa s, respectively. The relaxation time for 500 ppm, 1000 ppm, and 2000 ppm was calculated to be 9.1 ms, 12.4 ms, and 19.5 ms, respectively. The inlet boundary condition was set as fully developed flow and the outlets were set as 0 Pa for pressure control. The geometric domain was meshed with locally encrypted meshing at the outlets (shown in Figure 1c).

2.2. Theoretical Background

The viscoelastic effect was generated from the polymers added into the solutions. In the bounded flow field, when the viscoelastic fluids were experiencing shear flow, the polymer chains in the solutions were disturbed by the fluids and stretched along the main flow direction, leading to the anisotropy of the stress distribution. When there were suspended microparticles in viscoelastic fluids, the positive stresses on two sides of the microparticles were different, pushing the microparticles to migrate laterally while flowing downstream. Specifically, the microparticles in viscoelastic solutions performed a migration focus, which was subject to the elastic force pointing in the direction of the lowest shear rate.
For quantitative analysis, the force analysis of a certain micro-element showed that the micro-element was subjected to two normal stress differences, the first normal stress difference N 1 ( = τ xx τ yy ) and the second normal stress difference N 2 ( = τ yy τ zz ) . Generally, the first normal stress difference is positive and decreases as the solution volume fraction increases, while the second normal stress difference is negative and increases along with the solution volume fraction. For PEO viscoelastic solutions, N 2 is much smaller than N 1 and can be reasonably neglected. It is mainly the presence of the first normal stress difference N 1 that determines the lateral elastic forces on microparticles. The elastic force F E can be expressed as [39]:
F E a 3 N 1 = a 3 ( τ xx τ yy )
The strength of fluid viscoelasticity can be measured by a dimensionless Weissenberg number ( W i ) [40]:
W i = λ t f = λ γ ˙ = λ 2 U m w = 2 λ Q h w 2
where λ is the relaxation time of the fluid, U m and t f are the averaged velocity and characteristic time of the channel flow, respectively. The characteristic time is approximately equal to the inverse of the averaged (characteristic) shear rate γ ˙ , which is 2 U m / w or 2 Q / h w 2 in a rectangular channel.
The magnitude of the inertial effect is usually characterized using the Reynolds number, which is a dimensionless number that can be used to characterize the flow of a fluid. The expression of Reynolds number is as follows [41]:
Re = ρ U m D h μ = 2 ρ Q μ ( w + h )
where ρ is the density of the fluid, U m is the flow rate of the fluid, D h is the hydraulic diameter of the flow channel, μ is the dynamic viscosity. D h of the flow channel is related to the type of flow channel cross-section; for a rectangular cross-section flow channel commonly used in microfluidic control, the value of the hydraulic diameter can be estimated as D h = 2 w h / ( w + h ) , in which w and h are the width and height of the flow channel cross-section, respectively.
The vector of wall lift force and shear-gradient lift force are combined to be the net inertial lift force, which can be expressed as [42]:
F L = ρ U m 2 a 4 D h 2 f L ( R c , x c )
The magnitude of the dimensionless lift coefficient f L ( R c , x c ) is related to the position of the microparticles and the local Reynolds number.
A drag force arises when an object moves through a fluid or when the fluid flows past an object, due to a velocity difference between the particle and the fluid, which can be estimated by [38,43]:
F D = 3 π μ a ( v f v p )
where v f and v p are the velocities of the fluid element and particles, respectively.

2.3. Governing Equations

Based on the assumption that the viscoelastic fluid is incompressible and a continuous medium, the continuity equation and the momentum conservation equation in the control equations of the viscoelastic fluids are as follows [43,44]:
ρ u = 0
ρ ( u ) u = [ p I + K + T e ]
λ T e m + exp [ λ ε μ p tr ( T e m ) ] T e m = 2 u p D
D = 1 2 ( u + ( u ) T )
T e m = ( u ) T e m u T e m T e m ( u ) T
where T e = m T e m is the viscoelastic component of the stress tensor, which can be described by a different constitutive model [45,46,47]. In this paper, the first normal stress difference plays an important role in the flow field distribution of viscoelastic fluids. When establishing the connection between the normal stress difference and the high Wiesenberger number, an additional parameter ε is added to the Oldroyd-B model to create the PPT (Phan-Thien/Tanner) model. This additional parameter ensures the convergence of the model [44]. u is the velocity vector of the fluid, p is the pressure of the flow field, I is the unit tensor, K = 2 u s D is the Newtonian component of the stress tensor, u s and u p are solvent viscosity and polymer viscosity, which sum to u , and the retardation factor β is defined as β = u s / ( u s + u p ) . D is the strain velocity tensor, ε is the rheological parameter of the PPT model.

2.4. Boundary Conditions

Regarding boundary conditions, a fully developed flow characterized by a parabolic velocity profile is prescribed at the inlet:
v y = 0 , v x   is   parabolic   profile .
For channel flow, there is symmetry in the flow at the centerline. This results in both the normal flow and tangential stress becoming zero due to the symmetry condition present at the center line:
v x / y = 0 , v y = 0 , σ x y = 0
As the velocity profile and the values of stress tensor components are critical in resolving the given issue, particularly when they remain constant at a considerable distance from the inlet, we have established the stress tensor at the channel’s inlet:
σ v x x * = σ v x y * = σ v y y * = 0
The model incorporates no slip conditions for velocity at the channel wall. Furthermore, the boundary condition for the developed flow at the three outlets is established through the use of a pressure outlet condition:
[ p I + K + T e ] n = - p 0 n , p 0 = 0
To simulate viscoelastic fluid, the stationary solution is initialized before solving the problem. The commercial software COMSOL Multiphysics 6.0 utilizes the finite element method (FEM) to effectively solve governing equations. The PPT (Phan-Thien/Tanner) model is employed to obtain convergence and stability for large Weissenberg numbers. A PARDISO solver is preferred to the MUMPS solver for simulation with relative and absolute tolerances of 10−5 and 10−6, respectively, and a time step of 0.0001 s is used. In order to achieve convergence of the solution and determine the number of iterations required, the residual variables are adjusted until their values fall below the specified tolerance levels.
The COMSOL Multiphysics 6.0 software was utilized to perform a flow simulation lasting roughly 90 s, resulting in the acquisition of a stationary study of the fluid field. The subsequent particle tracking simulation utilized the previously determined steady-state solution as the initial condition to determine the particle migration trajectory. Each working condition was calculated for approximately 5 min.

2.5. Grid Independence

To assess the mesh independence of the numerical solution, four types of meshes—namely, regular, refined, more refined, and hyperfine—were employed. The number of grid cells was 13,513, 21,262, 32,288, and 135,078, with increasing density. The geometric meshes of the four different grids are shown in Figure 2a. The magnitudes and profiles of axial velocities near the inlet were compared at 500 ppm with an inlet flow rate of 30 μL/min for different mesh divisions, and the results are shown in Figure 2b. The axial velocity distribution results exhibit minimal discrepancies across various grids. To optimize the computational efficacy, this study opted for the refined grid, which boasts adequate precision.

3. Results and Discussion

3.1. Flow Field Analysis

The stationary flow field in the microchannel with a 500 ppm PEO solution was calculated under different flow rates from 10 μL/min to 80 μL/min. According to Equations (2) and (3), W i and R e will increase with the flow rates. Figure 3a shows a streamline distribution of the viscoelastic fluids flowing through the outlets under different flow rates. It can be observed that as the flow velocity increased, the flow streamlines gradually converged towards Outlet Ⅰ. This may result from the effects of additional cavities, and the increase in the flow rate magnified this differential. As the density of streamlines is proportional to the magnitude of the velocity at that point, this phenomenon indicates that the microparticles would be more likely to be dragged into Outlet Ⅰ. Under different flow conditions, the inlet velocity had a consistent variation pattern, where the maximum local velocity occurred at the center of the microfluidic channel (shown in Figure 3c). When the inlet velocity was 40 μL/min, the maximum velocity at the center of the microfluidic was 0.17 m/s, while it reached 0.35 m/s when the inlet velocity was 80 μL/min. The parameter ε serves as a restriction on the possible elongational viscosity values, meaning that an increase in this parameter results in a decrease in elongational viscosity while concurrently introducing elongational and shear thinning into the fluid model. It is evident that the PTT model can be simplified to the Oldroyd-B model by setting ε to zero. Figure 3c,d suggests that elevating the Weissenberg number and extensibility parameter ε results in a higher velocity gradient at the inlet. This outcome appears to be primarily due to the amplification of the shear-thinning characteristics inherent in PTT flow. The impact of the extensibility parameter ε on the flow field of this microfluidic channel is clearly less when juxtaposed with the influence of the Wiesenberger number. Consequently, for the purpose of this paper, the extensibility parameter ε was maintained at a constant value of 0.5.
Microparticles with the diameter of 4.8 µm were released in a viscoelastic fluid at 30 µL min−1 ( R e = 0.9, W i = 22.92). As a comparison, Newtonian flow under the same viscosity parameters was also simulated. The Newtonian flow presented a similar velocity distribution to the viscoelastic flow in the development phase, where the boundary layer of the microfluid developed sufficiently towards the center of the pipe. A fixed flow profile symmetrical to the center of the pipe was formed and the shape of the longitudinal component of the velocity vector approached a parabolic form, which corresponded to the flow of a Newtonian fluid. In the viscoelastic flow, the microparticles were deflected to the other side of the microchannel cavity and were collected from Outlet I, while in the Newtonian flow the microparticles were also deflected to some extent, but towards a lower level and were collected from Outlet II, as shown in Figure 4a. To explore the elastic effect of particle focusing, the microparticles were released at a fixed position, which was the center of the microchannel Inlet. The microparticles in the viscoelastic and Newtonian flows gradually migrated from the initial position to the opposite direction of the cavity. Figure 4b shows the particle trajectories within the two flow conditions over time. When the particle leaves the straight channel and enters the cavity, the lift force and drag force act together to move the particle toward the tip of the right triangle cavity. When the particle enters the cavity, the lift force and drag force change direction and, together with the elastic force, cause the particle to migrate in the direction away from the cavity. It is found that the microparticles in the viscoelastic flow exhibited a faster lateral migration rate than the Newtonian flow when moving into the same cavity under the same flow rate condition.
In the context of Newtonian flow, the particle undergoes both migration and rotation coinciding with the shear flow. The lateral particle migration is attributed to the inertial lift force, which encompasses both the wall-induced inertial lift and the shear-induced inertial lift [36]. While the wall-induced inertial lift functions to deflect the particle from the channel wall, the shear-induced inertial lift serves to push the particle towards the channel centerline. Through the balance between the wall-induced inertial lift and the shear-induced inertial lift, the particle is able to occupy a position somewhere between the centerline and the wall of the channel. In the viscoelastic flow, our numerical simulations reveal a notable occurrence: the stretching of particles within this type of flow generates an elastic force that prompts lateral migration towards the side that is distanced from the cavity. As a result, particles exhibit a considerably higher rate of lateral migration toward Outlet Ⅰ in viscoelastic flow when compared to Newtonian flow. The results indicated that the elastic effect of viscoelastic flow plays an important role in the lateral migration of microparticles.

3.2. Validation and Discussion of Microparticle Manipulation

Based on the obtained flow field calculations, a particle-tracing simulation module was integrated to predict the migration trajectory of microparticles of 4.8 μm throughout the microchannel. To calibrate and validate the proposed numerical model, the effects of PEO concentration and the flow rate on the focusing of 4.8 μm particles were investigated experimentally. The simulation results agreed well with the experimental results, as shown in Figure 5a. Figure 5b presents the various cases examined, where the PEO concentrations were 500 ppm, 1000 ppm, and 2000 ppm, and the flow rates were from 20 μL/min to 70 μL/min, respectively.
It can be derived from the PPT (Phan-Thien/Tanner) viscoelastic constitutive equation that the elastic effect of viscoelastic fluid in this model is determined by Wi, β , and ε . In this study, ε is a constant, and the changes in PEO solution concentration and the flow rate will result in the changes in W i and β . In a 500 ppm PEO solution at a flow rate of 30 μL/min ( β = 0.51, W i = 22.92), it was demonstrated in the calculated microparticle trajectories that the microparticles gradually migrate towards the opposite direction of the cavity when passing through the microchannel. When the flow rate increased to 40 μL/min ( β = 0.51, W i = 30.56), all the particles were tightly focused and flowed out from Outlet Ⅰ due to the balance of drag force and elastic force. Similar phenomena were observed in 1000 ppm and 2000 ppm PEO solutions, and the migration performance in 500 ppm PEO solution was better than that in 1000 ppm and 2000 ppm PEO solutions at the same flow rate. From Figure 6a, it can be seen that the flow rate is proportional to the particle lateral migration rate for the same PEO solution concentration. The PEO solution concentration is inversely proportional to the particle lateral migration rate for the same flow rate conditions. The optimal flow rate for 4.8 μm particle migration in 500 ppm ( β = 0.51, 30.56 ≤ W i ≤ 45.84) PEO solution was 40–60 μL/min.
To explore the physical mechanisms, Figure 6b shows a plot of the total force, elastic force, and drag force on the microparticles in the 500 ppm PEO solution versus the position of the triangle cavity at the flow rate of 50 µL/min. The force to the left of the green line is negative and points to the side of the cavity, while the force to the right of the green line is positive and points to the opposite side of the cavity, causing the particles to migrate towards Outlet Ⅰ. It can be seen that most of the total forces are located to the right of the green line, which is the root reason why the microparticles migrated laterally after passing through the cavities. The drag force on the particle changes significantly with the position, while the value of the elastic force is relatively stable and causes the particle to migrate towards the Outlet Ⅰ side.
The results suggest that the elastic effect of viscoelastic flow plays an important role in the lateral migration of particles. Figure 6c compares the elastic forces on the particles located in the middle position of the cavity in 500 ppm, 1000 ppm, and 2000 ppm PEO solutions at flow rates from 20 μL/min to 70 μL/min. It can be seen that with the increase in the flow velocity, the elastic force on the particles increases, while at the same flow velocity, the increase in PEO concentration will reduce the elastic force.
As per the formulation in Equation (1), the elastic force exhibits a correlation with the first normal stress difference N 1 . In the context of planar Poiseuille flow, wherein the velocity profile is curved and there exists a gradient in shear rate, the force attributable to the first normal stress difference N 1 exerts an influence in directing particles towards the central axis where the local shear rate is known to be the lowest [48]. The presented data in Figure 7a illustrate the distribution cloud map of the first normal stress difference N 1 observed in a 500 ppm PEO solution. It is apparent that the lowest shear rate is present at the cavity’s center. As the inlet flow rate escalates from 30 μL/min to 60 μL/min, the lowest shear rate progressively shifts away from the cavity towards Outlet I, and in tandem, the disparity between the maximum and minimum values of the first normal stress difference expands with the augmentation in flow velocity, resulting in a consequential surge in elastic force. The route of particle migration is depicted by the blue line. Meanwhile, this driving force depends strongly on particle size. As the particle diameter increases, a corresponding increase in elastic force contributes to a faster lateral migration rate, wherein the migration trajectories of particles measuring 3 μm, 4.8 μm, and 10 μm in diameter are determined using a calculation model in a 500 ppm PEO solution with an inlet flow rate of 30 μL/min. The calculation results align with the aforementioned rationale, as evidenced in Figure 7b.

4. Conclusions and Study Limitations

In this study, we established a 2D numerical model to simulate the microparticle trajectory in viscoelastic fluids and Newtonian fluids throughout a microchannel with additional triangular cavity arrays. The numerical model was calibrated and validated by the experimental results which showed good agreement. The mechanism of microparticle lateral migration was analyzed and discussed quantitatively based on the numerical results. It was observed that flow rate is proportional to particle lateral migration rate while PEO solution concentration is inversely proportional to particle lateral migration rate. This was because the elastic effect plays an important role in the lateral migration of particles, with an increase in flow velocity, PEO concentration, and elastic force on particles’ changes. It was observed that the first normal stress difference at peak stenosis increases with varying velocity of fluid for steady flow. The center of the cavity exhibits the lowest shear rate, and increasing the inlet flow rate causes the location of the lowest shear rate to shift towards Outlet I. As the flow velocity increases, the disparity between the values of the first normal stress difference widens, leading to a surge in elastic force. Based on the proposed numerical simulation methods and the findings of this study, one can undertake the optimization of the structural design of the microfluidic channels, identification of the optimal flow rate, and determination of the fluid physical characteristics such as viscoelasticity.
However, these findings need to be interpreted cautiously, and this study has several limitations. It is possible for the high Weissenberg number flow of a viscoelastic fluid in a contraction–expansion geometry to induce elastic turbulence; this induced turbulence can be further investigated, along with its destabilizing factors and characteristic captures, as well as work on the exploration of numerical stability. In addition, the viscoelastic flow behavior may change because each intrinsic constitutive model describes a particular fluid microstructure. As a result, by contrasting the variations in flow characterization, numerical simulations may be carried out for various intrinsic constitutive models to examine the impact of polymer microstructure on the flow and related rheological behavior.

Author Contributions

Conceptualization, T.W. and D.Y.; methodology, T.W.; software, B.Z.; validation, D.Y.; formal analysis, T.W. and W.W.; investigation, T.W.; resources, D.Y.; data curation, T.W., W.W., and D.Y.; writing—original draft preparation, T.W. and D.Y.; writing—review and editing, W.W. and B.Z.; visualization, T.W.; supervision, B.Z. and D.Y.; project administration, D.Y. All authors have read and agreed to the published version of the manuscript.

Funding

This work was supported by the National Natural Science Foundation of China (No. 52079122).

Data Availability Statement

Data will be made available on request.

Conflicts of Interest

The authors declare no conflict of interest.

Nomenclature

c [ppm]concentration
N1 [Pa]First normal stress difference
N2 [Pa]Second normal stress difference
WiWeissenberg number
Um [m/s] Averaged velocity
tf [s] Characteristic time
Q [m3]Average flow rate
ReReynolds Number
Dh [m]Hydraulic diameter
u [m/s]Velocity vector
p [Pa] Pressure
TeViscoelastic component of stress tensor
IUnit tensor
KNewtonian component of stress tensor
DStrain velocity tensor
λ [s] The relaxation time of fluid
τ [Pa] Stress tensor
γ [1/s]Shear rate
ρ [kg/m3]Density
μ [kg/ms] Dynamic viscosity
μs [kg/ms] Solvent viscosity
μp [kg/ms] Polymer viscosity
βRetardation factor
εRheological parameter of the PPT model
FE [N]Elastic force
FD [N]Lift force
FL [N]Drag force
w [μm]Width
h [μm]Height
a [μm]Particle diameter

References

  1. Hasani-Sadrabadi, M.M.; Taranejoo, S.; Dashtimoghadam, E.; Bahlakeh, G.; Majedi, F.S.; VanDersarl, J.J.; Janmaleki, M.; Sharifi, F.; Bertsch, A.; Hourigan, K.; et al. Microfluidic Manipulation of Core/Shell Nanoparticles for Oral Delivery of Chemotherapeutics: A New Treatment Approach for Colorectal Cancer. Adv. Mater. 2016, 28, 4134–4141. [Google Scholar] [CrossRef]
  2. Mohanty, S. Optically-actuated translational and rotational motion at the microscale for microfluidic manipulation and characterization. Lab Chip 2012, 12, 3624–3636. [Google Scholar] [CrossRef] [PubMed]
  3. Chen, P.; Li, S.; Guo, Y.; Zeng, X.; Liu, B.-F. A review on microfluidics manipulation of the extracellular chemical microenvironment and its emerging application to cell analysis. Anal. Chim. Acta 2020, 1125, 94–113. [Google Scholar] [CrossRef]
  4. Li, B.-L.; Li, D.-R.; Chen, J.-H.; Liu, Z.-Y.; Wang, G.-H.; Zhang, X.-P.; Xu, F.; Lu, Y.-Q. Hollow core micro-fiber for optical wave guiding and microfluidic manipulation. Sens. Actuators B Chem. 2018, 262, 953–957. [Google Scholar] [CrossRef]
  5. Sajeesh, P.; Sen, A.K. Particle separation and sorting in microfluidic devices: A review. Microfluid. Nanofluidics 2014, 17, 1–52. [Google Scholar] [CrossRef]
  6. Yeo, L.Y.; Chang, H.-C.; Chan, P.P.Y.; Friend, J.R. Microfluidic Devices for Bioapplications. Small 2011, 7, 12–48. [Google Scholar] [CrossRef] [PubMed]
  7. Kwizera, E.A.; Sun, M.; White, A.M.; Li, J.; He, X. Methods of Generating Dielectrophoretic Force for Microfluidic Manipulation of Bioparticles. ACS Biomater. Sci. Eng. 2021, 7, 2043–2063. [Google Scholar] [CrossRef]
  8. Wang, W.; Liu, Q.; Tanasijevic, I.; Reynolds, M.F.; Cortese, A.J.; Miskin, M.Z.; Cao, M.C.; Muller, D.A.; Molnar, A.C.; Lauga, E.; et al. Cilia metasurfaces for electronically programmable microfluidic manipulation. Nature 2022, 605, 681–686. [Google Scholar] [CrossRef]
  9. Zhao, W.; Cheng, R.; Miller, J.R.; Mao, L. Label-Free Microfluidic Manipulation of Particles and Cells in Magnetic Liquids. Adv. Funct. Mater. 2016, 26, 3916–3932. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  10. Zhang, Y.; Nguyen, N.-T. Magnetic digital microfluidics—A review. Lab Chip 2017, 17, 994–1008. [Google Scholar] [CrossRef] [Green Version]
  11. Wilson, R.; Reboud, J.; Bourquin, Y.; Neale, S.L.; Zhang, Y.; Cooper, J.M. Phononic crystal structures for acoustically driven microfluidic manipulations. Lab Chip 2011, 11, 323–328. [Google Scholar] [CrossRef] [PubMed]
  12. Zhang, P.; Bachman, H.; Ozcelik, A.; Huang, T.J. Acoustic Microfluidics. Annu. Rev. Anal. Chem. 2020, 13, 17–43. [Google Scholar] [CrossRef] [PubMed]
  13. Zheng, W.; Xie, R.; Liang, X.; Liang, Q. Fabrication of Biomaterials and Biostructures Based On Microfluidic Manipulation. Small 2022, 18, 2105867. [Google Scholar] [CrossRef] [PubMed]
  14. Reboud, J.; Bourquin, Y.; Wilson, R.; Pall, G.S.; Jiwaji, M.; Pitt, A.R.; Graham, A.; Waters, A.P.; Cooper, J.M. Shaping acoustic fields as a toolset for microfluidic manipulations in diagnostic technologies. Proc. Natl. Acad. Sci. USA 2012, 109, 15162–15167. [Google Scholar] [CrossRef] [Green Version]
  15. Zhu, P.; Wang, L. Passive and active droplet generation with microfluidics: A review. Lab Chip 2017, 17, 34–75. [Google Scholar] [CrossRef] [PubMed]
  16. Li, S.; Zhang, R.; Zhang, G.; Shuai, L.; Chang, W.; Hu, X.; Zou, M.; Zhou, X.; An, B.; Qian, D.; et al. Microfluidic manipulation by spiral hollow-fibre actuators. Nat. Commun. 2022, 13, 1331. [Google Scholar] [CrossRef] [PubMed]
  17. Guan, W.-S.; Huang, H.-X.; Chen, A.-F. Tuning 3D topography on biomimetic surface for efficient self-cleaning and microfluidic manipulation. J. Micromech. Microeng. 2015, 25, 035001. [Google Scholar] [CrossRef]
  18. Boran, Z.; Fan, Y.; Wenshuai, W.; Wuyi, W.; Wenhan, Z.; Qianbin, Z. Investigation of particle manipulation mechanism and size sorting strategy in a double-layered microchannel. Lab Chip 2022, 22, 4556–4573. [Google Scholar] [CrossRef]
  19. Lee, C.-Y.; Wang, W.-T.; Liu, C.-C.; Fu, L.-M. Passive mixers in microfluidic systems: A review. Chem. Eng. J. 2016, 288, 146–160. [Google Scholar] [CrossRef]
  20. Tofteberg, T.; Skolimowski, M.; Andreassen, E.; Geschke, O. A novel passive micromixer: Lamination in a planar channel system. Microfluid. Nanofluid. 2010, 8, 209–215. [Google Scholar] [CrossRef] [Green Version]
  21. Roudgar, M.; Brunazzi, E.; Galletti, C.; Mauri, R. Numerical Study of Split T-Micromixers. Chem. Eng. Technol. 2012, 35, 1291–1299. [Google Scholar] [CrossRef]
  22. Li, J.; Xia, G.; Li, Y. Numerical and experimental analyses of planar asymmetric split-and-recombine micromixer with dislocation sub-channels. J. Chem. Technol. Biotechnol. 2013, 88, 1757–1765. [Google Scholar] [CrossRef]
  23. Preira, P.; Grandné, V.; Forel, J.M.; Gabriele, S.; Camara, M.; Theodoly, O. Passive circulating cell sorting by deformability using a microfluidic gradual filter. Lab Chip 2013, 13, 161–170. [Google Scholar] [CrossRef] [PubMed]
  24. Tan, Y.-C.; Fisher, J.S.; Lee, A.I.; Cristini, V.; Lee, A.P. Design of microfluidic channel geometries for the control of droplet volume, chemical concentration, and sorting. Lab Chip 2004, 4, 292–298. [Google Scholar] [CrossRef] [PubMed]
  25. Yuan, D.; Zhao, Q.; Yan, S.; Tang, S.-Y.; Alici, G.; Zhang, J.; Li, W. Recent progress of particle migration in viscoelastic fluids. Lab Chip 2018, 18, 551–567. [Google Scholar] [CrossRef] [Green Version]
  26. Yuan, D.; Zhao, Q.; Yan, S.; Tang, S.-Y.; Zhang, Y.; Yun, G.; Nguyen, N.-T.; Zhang, J.; Li, M.; Li, W. Sheathless separation of microalgae from bacteria using a simple straight channel based on viscoelastic microfluidics. Lab Chip 2019, 19, 2811–2821. [Google Scholar] [CrossRef] [PubMed]
  27. Yuan, D.; Yan, S.; Zhang, J.; Guijt, R.M.; Zhao, Q.; Li, W. Sheathless Separation of Cyanobacterial Anabaena by Shape Using Viscoelastic Microfluidics. Anal. Chem. 2021, 93, 12648–12654. [Google Scholar] [CrossRef] [PubMed]
  28. Zhang, T.; Liu, H.; Okano, K.; Tang, T.; Inoue, K.; Yamazaki, Y.; Kamikubo, H.; Cain, A.K.; Tanaka, Y.; Inglis, D.W.; et al. Shape-based separation of drug-treated Escherichia coli using viscoelastic microfluidics. Lab Chip 2022, 22, 2801–2809. [Google Scholar] [CrossRef]
  29. Feng, H.; Jafek, A.R.; Wang, B.; Brady, H.; Magda, J.J.; Gale, B.K. Viscoelastic Particle Focusing and Separation in a Spiral Channel. Micromachines 2022, 13, 361. [Google Scholar] [CrossRef]
  30. D’Avino, G.; Maffettone, P.L. Particle dynamics in viscoelastic liquids. J. Non-Newton. Fluid Mech. 2015, 215, 80–104. [Google Scholar] [CrossRef]
  31. Huang, P.; Joseph, D. Effects of shear thinning on migration of neutrally buoyant particles in pressure driven flow of Newtonian and viscoelastic fluids. J. Non-Newton. Fluid Mech. 2000, 90, 159–185. [Google Scholar] [CrossRef]
  32. Seo, K.W.; Byeon, H.J.; Huh, H.K.; Lee, S.J. Particle migration and single-line particle focusing in microscale pipe flow of viscoelastic fluids. RSC Adv. 2014, 4, 3512–3520. [Google Scholar] [CrossRef] [Green Version]
  33. Li, G.; McKinley, G.H.; Ardekani, A.M. Dynamics of particle migration in channel flow of viscoelastic fluids. J. Fluid Mech. 2014, 785, 486–505. [Google Scholar] [CrossRef] [Green Version]
  34. Yang, S.; Kim, J.Y.; Lee, S.J.; Lee, S.S.; Kim, J. Sheathless elasto-inertial particle focusing and continuous separation in a straight rectangular microchannel. Lab Chip 2011, 11, 266–273. [Google Scholar] [CrossRef] [PubMed]
  35. Del Giudice, F.; D’Avino, G.; Greco, F. Effect of fluid rheology on particle migration in a square-shaped microchannel. Microfluid. Nanofluid. 2015, 19, 95–104. [Google Scholar] [CrossRef] [Green Version]
  36. Jiang, D.; Ni, C.; Tang, W.; Xiang, N. Numerical simulation of elasto-inertial focusing of particles in straight microchannels. J. Phys. D Appl. Phys. 2020, 54, 065401. [Google Scholar] [CrossRef]
  37. Naderi, M.M.; Barilla, L.; Zhou, J.; Papautsky, I.; Peng, Z. Elasto-Inertial Focusing Mechanisms of Particles in Shear-Thinning Viscoelastic Fluid in Rectangular Microchannels. Micromachines 2022, 13, 2131. [Google Scholar] [CrossRef] [PubMed]
  38. Yuan, D.; Zhang, J.; Sluyter, R.; Zhao, Q.; Yan, S.; Alici, G.; Li, W. Continuous plasma extraction under viscoelastic fluid in a straight channel with asymmetrical expansion–contraction cavity arrays. Lab Chip 2016, 16, 3919–3928. [Google Scholar] [CrossRef]
  39. Lim, H.; Back, S.M.; Hwang, M.H.; Lee, D.H.; Choi, H.; Nam, J. Sheathless High-Throughput Circulating Tumor Cell Separation Using Viscoelastic non-Newtonian Fluid. Micromachines 2019, 10, 462. [Google Scholar] [CrossRef] [Green Version]
  40. Di Carlo, D. Inertial microfluidics. Lab Chip 2009, 9, 3038–3046. [Google Scholar] [CrossRef]
  41. Asmolov, E.S. The inertial lift on a spherical particle in a plane Poiseuille flow at large channel Reynolds number. J. Fluid Mech. 1999, 381, 63–87. [Google Scholar] [CrossRef]
  42. Zhang, J.; Yan, S.; Yuan, D.; Alici, G.; Nguyen, N.-T.; Ebrahimi Warkiani, M.; Li, W. Fundamentals and applications of inertial microfluidics: A review. Lab Chip 2016, 16, 10–34. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  43. Tang, D.; Marchesini, F.H.; D’hooge, D.R.; Cardon, L. Isothermal flow of neat polypropylene through a slit die and its die swell: Bridging experiments and 3D numerical simulations. J. Non-Newton. Fluid Mech. 2019, 266, 33–45. [Google Scholar] [CrossRef]
  44. Xue, S.C.; Phan-Thien, N.; Tanner, R.I. Numerical study of secondary flows of viscoelastic fluid in straight pipes by an implicit finite volume method. J. Non-Newton. Fluid Mech. 1995, 59, 191–213. [Google Scholar] [CrossRef]
  45. Marchal, J.M.; Crochet, M.J. Hermitian finite elements for calculating viscoelastic flow. J. Non-Newton. Fluid Mech. 1986, 20, 187–207. [Google Scholar] [CrossRef]
  46. Brooks, A.N.; Hughes, T.J.R. Streamline upwind/Petrov-Galerkin formulations for convection dominated flows with particular emphasis on the incompressible Navier-Stokes equations. Comput. Methods Appl. Mech. Eng. 1982, 32, 199–259. [Google Scholar] [CrossRef]
  47. Marchal, J.M.; Crochet, M.J. A new mixed finite element for calculating viscoelastic flow. J. Non-Newton. Fluid Mech. 1987, 26, 77–114. [Google Scholar] [CrossRef]
  48. Zhou, J.; Papautsky, I. Viscoelastic microfluidics: Progress and challenges. Microsyst. Nanoeng. 2020, 6, 113. [Google Scholar] [CrossRef]
Figure 1. The schematics of the microchannel with asymmetrical expansion–contraction cavity arrays: (a) flow velocity distribution at the inlet; (b) dimensions of asymmetrical additional cavity arrays; (c) local encrypted meshing at the Outlets.
Figure 1. The schematics of the microchannel with asymmetrical expansion–contraction cavity arrays: (a) flow velocity distribution at the inlet; (b) dimensions of asymmetrical additional cavity arrays; (c) local encrypted meshing at the Outlets.
Micromachines 14 00915 g001
Figure 2. (a) Reticulation of geometry for four various meshes; (b) axial velocity profile near the inlet at 500 ppm with an inlet flow rate of 30 μL/min.
Figure 2. (a) Reticulation of geometry for four various meshes; (b) axial velocity profile near the inlet at 500 ppm with an inlet flow rate of 30 μL/min.
Micromachines 14 00915 g002
Figure 3. (a) Streamline distribution at different flow rates; (b) the simulated flow field at the outlet under different velocities; (c) comparison of flow velocity distribution at the inlet of the microchannel in viscoelastic fluids at different Weissenberg numbers ( W i ); (d) comparison of flow velocity distribution at the inlet of the microchannel for different ε ( R e = 0.9, W i = 22.92).
Figure 3. (a) Streamline distribution at different flow rates; (b) the simulated flow field at the outlet under different velocities; (c) comparison of flow velocity distribution at the inlet of the microchannel in viscoelastic fluids at different Weissenberg numbers ( W i ); (d) comparison of flow velocity distribution at the inlet of the microchannel for different ε ( R e = 0.9, W i = 22.92).
Micromachines 14 00915 g003
Figure 4. Simulation results of microparticle manipulation in viscoelastic fluid and Newtonian fluid. (a) Calculated results of microparticle migration in Newtonian and viscoelastic flows (the coordinate origin (0, 0) of the microchannel model is located at the upper vertex of the inlet). (b) Local migration trajectories of microparticles in Newtonian and viscoelastic flows.
Figure 4. Simulation results of microparticle manipulation in viscoelastic fluid and Newtonian fluid. (a) Calculated results of microparticle migration in Newtonian and viscoelastic flows (the coordinate origin (0, 0) of the microchannel model is located at the upper vertex of the inlet). (b) Local migration trajectories of microparticles in Newtonian and viscoelastic flows.
Micromachines 14 00915 g004
Figure 5. (a) Experimental results under microscopes and numerical simulation results of the 4.8 μm fluorescent microparticle distribution in the outlet section (1000 ppm PEO solution): (i) flow rate = 30 μL/min ( W i = 31 ), (ii) flow rate = 40 μL/min ( W i = 41.33 ); (b) the distribution of 4.8 µm microparticles in the outlet region under flow rates from 20 µL/min to 70 µL/min (500 ppm PEO solutions, 1000 ppm, and 2000 ppm).
Figure 5. (a) Experimental results under microscopes and numerical simulation results of the 4.8 μm fluorescent microparticle distribution in the outlet section (1000 ppm PEO solution): (i) flow rate = 30 μL/min ( W i = 31 ), (ii) flow rate = 40 μL/min ( W i = 41.33 ); (b) the distribution of 4.8 µm microparticles in the outlet region under flow rates from 20 µL/min to 70 µL/min (500 ppm PEO solutions, 1000 ppm, and 2000 ppm).
Micromachines 14 00915 g005
Figure 6. (a) Comparison of particle migration trajectories at different flow rates and PEO solution concentrations: (i) 500 ppm PEO solution with flow rates of 30 μL/min, 50 μL/min, and 70 μL/min, (ii) particle trajectories in 500 ppm, 1000 ppm, and 2000 ppm PEO solutions at a flow rate of 50 μL/min; (b) the total force, elastic force, and traction force of the particles in the 500 ppm PEO solution versus the position of the triangle cavity at a flow rate of 50 µL/min; (c) the elastic forces on the particles located in the middle position of the cavity in 500 ppm, 1000 ppm, and 2000 ppm PEO solutions at flow rates from 20 μL/min to 70 μL/min.
Figure 6. (a) Comparison of particle migration trajectories at different flow rates and PEO solution concentrations: (i) 500 ppm PEO solution with flow rates of 30 μL/min, 50 μL/min, and 70 μL/min, (ii) particle trajectories in 500 ppm, 1000 ppm, and 2000 ppm PEO solutions at a flow rate of 50 μL/min; (b) the total force, elastic force, and traction force of the particles in the 500 ppm PEO solution versus the position of the triangle cavity at a flow rate of 50 µL/min; (c) the elastic forces on the particles located in the middle position of the cavity in 500 ppm, 1000 ppm, and 2000 ppm PEO solutions at flow rates from 20 μL/min to 70 μL/min.
Micromachines 14 00915 g006
Figure 7. (a) The distribution cloud map of the first normal stress difference observed in a 500 ppm PEO solution; (b) comparison of particle migration trajectories at different particle diameters: 500 ppm PEO solution with flow rates of 30 μL/min, particle diameters are 3 μm, 4.8 μm, and 10 μm.
Figure 7. (a) The distribution cloud map of the first normal stress difference observed in a 500 ppm PEO solution; (b) comparison of particle migration trajectories at different particle diameters: 500 ppm PEO solution with flow rates of 30 μL/min, particle diameters are 3 μm, 4.8 μm, and 10 μm.
Micromachines 14 00915 g007
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Wang, T.; Yuan, D.; Wan, W.; Zhang, B. Numerical Study of Viscoelastic Microfluidic Particle Manipulation in a Microchannel with Asymmetrical Expansions. Micromachines 2023, 14, 915. https://doi.org/10.3390/mi14050915

AMA Style

Wang T, Yuan D, Wan W, Zhang B. Numerical Study of Viscoelastic Microfluidic Particle Manipulation in a Microchannel with Asymmetrical Expansions. Micromachines. 2023; 14(5):915. https://doi.org/10.3390/mi14050915

Chicago/Turabian Style

Wang, Tiao, Dan Yuan, Wuyi Wan, and Boran Zhang. 2023. "Numerical Study of Viscoelastic Microfluidic Particle Manipulation in a Microchannel with Asymmetrical Expansions" Micromachines 14, no. 5: 915. https://doi.org/10.3390/mi14050915

APA Style

Wang, T., Yuan, D., Wan, W., & Zhang, B. (2023). Numerical Study of Viscoelastic Microfluidic Particle Manipulation in a Microchannel with Asymmetrical Expansions. Micromachines, 14(5), 915. https://doi.org/10.3390/mi14050915

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop