Confinement effect on the viscoelastic particle ordering in microfluidic flows: Numerical simulations and experiments

The formation of self-assembled organized structures in microfluidic flows, also called microfluidic crystals, is not only fascinating under a fundamental point of view, but it is also linked to a plethora of applications across a diversity of fields, including biomedical engineering and materials science.^1,2 For instance, crystals made of highly confined droplets (droplet crystals) flowing in a microchannel develop vibrational modes reminiscent of the “phonons” in solid-state physics.³ Such behavior is generated by the long-range hydrodynamic interactions occurring in confined flows along the equally spaced droplets arranged in a one-dimensional (1D) array.⁴ In addition to 1D and 2D planar configurations, droplet crystals are formed in a 3D configuration too by designing an appropriate microfluidic device.^5,6 The understanding of the droplet crystals formation dynamics has also a great impact on the design of microfluidic systems for biomedical applications such as polymerase chain reaction (PCR)⁷ and determination of the minimum inhibitory concentration (MIC) for antibiotics⁸ and is relevant for the synthesis of novel complex materials.^9,10

Microfluidic crystals made of rigid particles have recently attracted interest for several applications related to encapsulation^11–13 and flow cytometry.¹⁴ Rigid particles form self-ordered structures, also called particle trains, resulting from the flow-mediated hydrodynamic interactions between consecutive particles.¹ The two main examples are inertially mediated and viscoelasticity-mediated hydrodynamic interactions responsible for the formation of inertial and viscoelastic particle trains, respectively.¹ Inertial particle trains were first observed by Segré and Silberberg¹⁵ for flow in a macroscopic pipe. Only many years later, Di Carlo et al.¹⁶ investigated the formation of inertial particle trains in microfluidic geometries. Since then, the mechanism driving the formation of trains of particles due to inertial effects has been extensively studied through experiments and numerical simulations;^17–27 however, a complete understanding of this phenomenon remains currently out of reach.

Similarly, the physical mechanism responsible for viscoelastic particle trains is still elusive, with only a few studies focusing on the subject. The potential existence of viscoelastic particle trains was first theorized by D'Avino et al.²⁸ They performed numerical simulations on pairs and triplets of particles, finding that particles experienced either attractive or repulsive forces depending on the flow conditions and the initial distance between consecutive particles. They suggested that the repulsion dynamics between consecutive particles, for a system of several particles, can lead to the formation of a particle train. The first experimental evidence of this phenomenon was reported in Del Giudice et al.,²⁹ where particle trains aligned on the channel centerline were observed when the suspending liquid presented elasticity and shear-thinning features (i.e., the viscosity of the fluid decreases when increasing the shear rate). In the same study, strings of attached particles were instead observed for a fluid with negligible shear-thinning. Liu et al.³⁰ proposed a microfluidic device to control the spacing between consecutive particles consisting of a first element aimed at breaking particle aggregates followed by a system to make the local concentration of particles uniform. More recently, particle trains were observed in different types of viscoelastic liquids³¹ and exploited to improve the encapsulation efficiency of particles above the Poisson stochastic value.¹³ Numerical simulations aiming at elucidating the dynamics of viscoelastic particle trains are also limited. D'Avino and Maffettone³² found that viscoelastic particle train formation in a cylindrical microchannel strongly depends on the initial particle distribution as well as on the particle volume fraction, while the volumetric flow rate has a small effect on the ordering efficacy. In agreement with experimental findings,^29,31 strings of particles in contact were observed at non-negligible volume fractions. Hu et al.³³ performed numerical simulations to evaluate the effect of several parameters on the formation of particle trains in power-law fluids, observing a rich variety of inter-related effects contributing to the particle train formation. Very recently, Hu et al.³⁴ presented a numerical study featuring the self-assembly dynamics of particle chains in a straight microchannel, observing that viscoelastic liquids with a negligible shear-thinning were responsible for the particle chaining, as also previously observed both numerically²⁸ and experimentally.²⁹ 

Although the aforementioned works clearly show the possibility to form self-ordered particle trains through fluid viscoelasticity, several outstanding points remain. A critical issue is the formation of strings of the particle in contacts that are detrimental for particle ordering. Furthermore, it is not clear the effect of the confinement ratio, defined as the ratio between the particle and channel size, on the self-ordering dynamics, being all the previous works carried out for one particle dimension. Understanding such effect is important for the design and optimization of microfluidic devices for a diverse set of applications, while also contributing to clarify the origin and impact of viscoelasticity-mediated hydrodynamic interactions on the particle self-ordering phenomenon.

In this work, we performed numerical simulations and experiments to elucidate the effect of the confinement ratio on the self-assembly dynamics of particle trains in viscoelastic liquids flowing in a straight microchannel. We first performed numerical simulations for different particle sizes and concentrations and evaluated the distributions of the interparticle spacing at several distances from the channel inlet. We then performed experiments using an aqueous solution of xanthan gum (XG) 0.2 wt. % for two different particle sizes, distances from the channel inlet and volumetric flow rates. We also experimentally demonstrated that a microfluidic device made of a series of elements with trapezoidal shape was capable of significantly reducing the amount of particle aggregates, thus leading to a more uniform particle train. Our results shed some light on the particle ordering phenomenon in viscoelastic liquids.

The methodology employed to simulate the multiparticle system follows the one proposed in D'Avino and Maffettone.³² The particles were assumed spherical with diameter d, suspended in a viscoelastic fluid, and aligned along the centerline of a cylindrical channel with diameter D. A system made of three particles was first investigated through direct numerical simulations, i.e., by solving the continuity and momentum balance equations for the fluid and the force balance for the particle assuming inertialess and force-free conditions. The liquid was modeled by the Giesekus constitutive equation.³⁵ A flow rate Q was imposed at the channel inlet, and periodic boundary conditions were applied between the inlet and outlet channel sections. No-slip conditions were set at the channel wall and at the particle surfaces, resulting in the rigid-body motion condition. An Arbitrary Lagrangian–Eulerian finite element method was used to solve the governing equations and handle the particle motion.

The results for the three-particle system were used to build a lookup table that can be interpolated to obtain the particle velocities for specific interparticle distances.³² The basic assumption for the multiparticle system simulation was that each particle hydrodynamically interacts with the previous and next one only. Hence, the train can be decomposed in a set of several systems made by three particles with velocities taken from the lookup table. The motion of each particle can then be computed by integrating the motion equation $dzi/dt=Vi(s1,i,s2,i)$⁠, where z_i and V_i are the axial coordinate and velocity of the $i−$ th particle, where V_i depends on the distances $s1,i$ and $s2,i$ between the particle and the previous and next one. More details on the mathematical model, the adopted numerical method, and the simulation procedure can be found elsewhere.³² 

An aqueous xanthan gum (XG) solution at 0.2 wt. % was prepared by dissolving Xanthomonas campestris (Sigma Aldrich, UK) in deionized water. XG solutions were previously employed in particle ordering experiments.^13,31 The solution was mixed using a magnetic stirrer for 24 h to allow full dissolution of the polymer. The rheological characterization was carried out using a stress-controlled rheometer (TA AR2000ex) with an acrylic cone (60 mm diameter, 1° angle) at a constant temperature $T=20 °$C. XG solution displayed a zero-shear viscosity plateau at shear rate values $γ̇<1 s−1$⁠, followed by a shear-thinning region at $γ̇>1 s−1$ [Fig. 1(a)], in line with previous experimental characterizations on the same material.^13,31,36 Small amplitude oscillatory shear (SAOS) was used to evaluate the fluid viscoelasticity, quantified via the storage modulus $G′$ and the loss modulus $G″$ [Fig. 1(b)]. From the SAOS experiments, we could not evaluate directly the longest relaxation time λ, because we did not observe the terminal region at low angular frequencies. Following the same procedure employed previously,^13,31,36 we identified the intersection point between the power-law fit and straight line representing the zero shear viscosity plateau to determine the critical shear rate $γ̇c$ for the onset of the shear-thinning and then evaluated the longest relaxation time as $λ=1/γ̇c≃1.25$ s.

For the microfluidic experiments, rigid polystyrene particles (Polysciences) having diameters of $d=20±2$ μm and $d=45±5$ μm were added to the XG solutions at concentrations of 0.2 wt. % and 0.3 wt. % for the d = 20 μm particles and 0.7 wt. % for the d = 45 μm particles. A vortex stirrer was used to disperse the particle uniformly in the XG suspending fluid.

We performed experiments on two different microfluidic devices. The first one consisted of a series of 16 trapezoidal elements with sections increasing along the flow direction as displayed in Fig. 2(a), and it was used to break particle aggregates formed at the device inlet resulting from the relatively large number of particles, similar to previous studies.³⁰ The second device was a combination of the first device followed by a straight 30 cm long glass capillary with a circular section employed to experimentally evaluate the evolution of the particle self-assembly along the channel axis [Fig. 2(b)].

A micromilling machine (Minitech CNC Mini-Mill) was employed to fabricate the microfluidic device molds on a rigid poly(methyl methacrylate) (PMMA) substrate with a thickness of 1.2 mm. First, a 2-mm-wide tip was used to mill all the substrate down to 300 μm to achieve a uniform surface. Then, a 1.1 mm tip was used to mill the channel onto the PMMA substrate, while a 100 μm tip was used to mill around the sharp edges of the trapezoidal elements to maintain an acute angle. The channel depth was kept constant to 100 μm for the device in Fig. 2(a) and 175 μm for the device in Fig. 2(b): this was due to the fact that the external diameter of the glass capillary was 175 μm. After fabrication, the devices were placed in an ultrasonic bath for 15 min to remove the excess material. Pressurized air was used to dry the PMMA substrates. A polydymethyl siloxane silicone elastomer and curing agent (SLYGARD) prepared at a ratio of 1:10 were mixed thoroughly in a falcon tube. The viscous solution was poured onto the PMMA mold. A desiccator was employed to remove air bubbles. The chamber was then placed in an oven (Fisher Brand) pre-heated to 72 °C for 1 h to allow the polydimethylsiloxane (PDMS) to cross-link. A sharp blade was used to cut the PDMS device, before placing the device in an ultrasonic bath for 10 min to remove any impurity. A plunger with a diameter of 1.2 mm was used to create the inlet and the outlet of the device. A double-sided tape (Adhesives Research, Inc.) was employed to bond the PDMS device to a glass slide. A round glass capillary (Vitrocom, CM Scientific) with an internal diameter of 100 μm and an external diameter equal to 175 μm was placed directly into the PDMS device, as schematized in Fig. 2(b). A viscous glue (Gorilla Inc.) was put onto the interface between the PDMS mold and the glass capillary to prevent fluid leaking from the PDMS. The glue was left to dry for 24 h.

An inverted microscope (Zeiss Axiovert 135) connected to a fast camera (Photron Mini UX50) was used to record videos of particles flowing in the two microfluidic devices described above. The flow was controlled using a pressure pump (Dolomite Microfluidics) connected to an external compressor (Balbi) for air supply. For all the experiments with the 20 μm particles, a pressure drop $ΔP=1000$ mbar was initially imposed to flow the suspension into the channel. Then, the pressure drop was lowered to $ΔP=$ 600 mbar, 500 mbar, 400 mbar, and 300 mbar to record the videos. For the 45 μm particle experiments, an initial pressure drop $ΔP=2000$ mbar was imposed to flow the suspension into the channel, before lowering it to $ΔP=$ 1000 mbar, 900 mbar, 800 mbar, and 700 mbar to record the videos. Flow rate values were evaluated using a flow sensor (Dolomite Microfluidics), with values reported in Table I. Particle tracking was carried out using a Matlab program employed before.³¹ The particles center-to-center distance, s, was divided by the particle diameter to obtain the normalized distance between particles $S*=s/d$ as done in previous works.^13,29,31 The frequency distribution for $S*$ data was computed using Origin Pro. The bin size for $S*$ was set equal to 1, corresponding to one particle diameter.

The effect of the confinement ratio $β=d/D$⁠, representing the ratio between the particle diameter d and channel diameter D, on particle ordering was first addressed by numerical simulations for a three-particle system schematized in Fig. 3(a). For a fixed value of the confinement ratio, several simulations by varying the two initial relative particle distances [s₁ and s₂ in Fig. 3(a)] were carried out. The simulations were run for a sufficiently long time such that the viscoelastic stresses fully develop. At the end of each simulation, the relative distances between the particles and the corresponding velocities were stored and used to build a lookup table for the multiparticle system simulations, similar to the approach employed before.³² All the simulations were done by fixing the following dimensionless parameters: De = 2, $ηr=0.1, α=0.2$⁠, where $De=4λQ/(πD3)$ is the ratio between the fluid relaxation time λ and the flow characteristic time, $ηr=ηs/ηp$ is the ratio between the (Newtonian) solvent and the polymer viscosity, and α is a constitutive parameter. The chosen value for α gives a shear-thinning viscosity in the shear flow. For this set of parameters, the dynamics between two particles is mainly repulsive regardless of the interparticle distance and the confinement ratio.³⁷ This repulsive pairwise interaction was observed in a wide range of flow rates, confinement ratios, and fluid rheology,³⁷ and it is the basic ingredient to explain the formation of stable trains.²⁹ Hence, we expect that the train dynamics obtained with these parameters is representative of a variety of fluids and flow conditions and should promote the formation of ordered structures.

The analysis of the three-particle system already revealed useful information on the effect of the confinement ratio on the particle dynamics. With reference to Fig. 3(a), the trailing, middle, and leading particles are labeled as “1,” “2,” and “3,” respectively, s₁ and s₂ are the distances between the particle centers, and $σ1=s1−d$ and $σ2=s2−d$ are the distances between the particle surfaces. The panels in Fig. 3(b) show the relative velocity between the middle and the trailing particles, $ΔV21=V2−V1$⁠, as a function of the distance σ₂ for fixed values of the distance σ₁, where both distances were made dimensionless by the channel diameter. (The dimensionless distances are denoted by a superscript “D,” i.e., $σ1D=σ1/D$ and $σ2D=σ2/D$⁠.) In each plot, five curves are displayed corresponding to different values of the confinement ratio $β=d/D$⁠. Due to the definition of σ, all the curves start from 0 corresponding to the middle and leading particles in contact. Depending on the two distances, the relative velocity can be negative or positive denoting that the middle and leading particles can approach or separate while traveling along the channel. The value of $σ2D$ such that the relative velocity changes sign was independent of the confinement ratio. Also note that, when the trailing particle was far from the middle one (see the bottom panel at $σ1D=3$⁠), the system reduced to two particles, where the dynamics was repulsive for any interparticle distance, in agreement with previous simulation results.³⁷ The confinement ratio, however, has a relevant quantitative effect on the relative particle motion. By reducing the confinement from $β=0.4$ to $β=0.2$⁠, the relative velocity decreased by approximately one order of magnitude. This is related to the weaker hydrodynamic interactions as the particle size is reduced. Hence, for small confinement ratios, a much longer microchannel is needed to observe a variation of the interparticle distances as compared to more confined geometries.

The curves at different confinement ratios can be scaled using $β4$ as scaling factor as visible in the panels of Fig. 3(c), where the relative velocities are multiplied by $(βref/βi)4$ with $βref$ a reference value and β_i the values of the confinement ratio corresponding to the different curves. We chose $βref=0.4$⁠, so all the curves collapsed on the purple one. The plots clearly indicate that the scaling works well as the particles were not very close to each other, i.e., for $σ2D>0.5$⁠. Indeed, for closer particles, the scaled curves at different confinement ratios deviated, and we were unable to find a simple scaling law in this region. A similar scaling was found as the distance $σ1D$ was varied and for the relative velocity between the leading and middle particles $ΔV32=V3−V2$ (the corresponding curves are reported in Figs. S1 and S2). In conclusion, the dynamics of a triplet of particles with different confinement ratios, but at the same relative distances between the particle surfaces $σD>0.5$⁠, follows the same dynamics. Hence, the evolution of the interparticle distances observed for a strong confined geometry will be observed for a small confinement ratio as well as at a larger distance from the channel inlet. We also tried to plot the relative velocities as a function of the distance between the particle surfaces made dimensionless by the particle diameter, i.e., $σd=σ/d$⁠. Note that, according to this kind of a dimensionalization, the train microstructure at two different confinement ratios is self-similar, i.e., it is simply rescaled. This was used in Del Giudice et al.²⁹ to simulate the train dynamics at $β=0.2$ from data at $β=0.4$⁠. When plotted in terms of σ^d, the relative velocity curves can be roughly scaled by using the same factor as before with a third power instead of four, as assumed in Del Giudice et al.²⁹ on the basis of results for a particle pair. However, the scaling was not as good as the one displayed in Fig. 3(c), suggesting that the relevant dimension is the distance between the particle surfaces relative to the channel size rather than to the particle diameter, as also found for inertial ordering.²² 

As previously mentioned, the results for the three-particle system were stored in a lookup table used to evaluate the dynamics of a particle train.³² Figure 4 reports the distributions of the interparticle distances for different values of the confinement ratio, linear concentration $ϕ$⁠, and several distances from the channel inlet L/D with L the channel length. The linear concentration is defined as the length of the channel centerline occupied by the particles divided by the total channel length, i.e., $ϕ=Nd/L$⁠, where N is the number of particles flowing in the channel with length L. All the simulations considered a uniform distribution in inflow (see the red bars in each panel). To compare the distributions at different confinement ratios and linear concentration, the following normalized distance is defined:³² 

where σ₀ is the average of the distances between the surfaces of two consecutive particles at the channel inlet. In this way, the initial distances were uniformly distributed between 0 and 1 regardless of the confinement ratio and the linear concentration.

By varying the confinement ratio, similar features were observed: (i) at low values of $ϕ$⁠, a preferential distance was visible, denoting the tendency of the train to self-assemble in an ordered structure; (ii) as $ϕ$ increased, the number of strings of attached particles increased (see the bars at the lowest value of $σnorm$⁠) in agreement with previous works;^29,31 and (iii) at the highest investigated linear concentration $ϕ=0.4$⁠, ordering was suppressed and the particles were arranged in strings made of several particles in contact. As expected, the arrangement of the microstructure was slower for low $β$–values. Specifically, according to the three-particle system results discussed above, the self-similar train microstructure at two different confinement ratio values should roughly follow the same dynamics scaled by a factor that depends on the confinement ratio to the third power. Considering that a variation of the confinement does not significantly alter the train average velocity (the deviations of the axial velocity of an isolated particle at $β=0.2$ and $β=0.4$ are less than 3%), we expect that similar interparticle distance distributions are observed further from the channel inlet as β is reduced. As previously remarked, this behavior is expected as the particles are sufficiently separated in inflow (⁠$σD>0.5$⁠), which more likely occurs at low linear concentrations. Indeed, the histograms of the interparticle distances for $β=0.2$ and $L/D=3000$ were similar to those at $β=0.3$ and $L/D=1000$ and at $β=0.4$ and $L/D=250$ (see Fig. S3). If the linear concentration is large, the probability to find close particles (⁠$σD<0.5$⁠) increases, and the scaling is lost. Indeed, the formation of particle strings is a fast phenomenon regardless of the confinement ratio as clearly visible from the rapid increase of the bars at low values of $σnorm$⁠, which is quantitatively similar to any $β$-value when $ϕ=0.3$ and $ϕ=0.4$⁠.

The numerical simulations highlighted three important aspects of the particle ordering phenomenon. First, for a fixed tube diameter, larger particles tend to self-order on shorter distances compared to smaller particles. Second, the ordering dynamics is strongly controlled by the distance between the surfaces of consecutive particles. Third, large particle concentration results in the formation of particle aggregates that may suppress the formation of an equally spaced microstructure. We now present the experimental results performed using an aqueous xanthan gum (XG) 0.2 wt. % solution, similar to previous ordering studies.^13,31 The aqueous XG 0.2 wt. % is elastic and shear-thinning (Fig. 1), two attributes required to investigate the validity of the numerical simulations results performed using the Giesekus model at De = 2, i.e., when the fluid experiences a shear-thinning behavior.

The initial distribution of the distances between consecutive particles entering the channel has a strong influence on the train dynamics affecting, in turn, the ordering efficiency.^31,32 Previous numerical³² and experimental^29,31 results demonstrated that the presence of particle aggregates, e.g., doublets and triplets of particles, is detrimental to the formation of equally spaced structures, thus reducing the number of particle aggregates in inflow is essential for the formation of a uniform particle train. In our previous work,³¹ we showed that limiting the number of connections between the entrance of the microfluidic device and the pump resulted in an overall reduction of particle strings. However, this remedy is not always applicable, and the best solution would actually be the design of a microfluidic device capable of breaking particle aggregates before they start to hydrodynamically interact and form self-ordered structures. We designed a microfluidic channel similar to the one proposed by Liu et al.,³⁰ made of a series of 16 elements with trapezoidal shape [Fig. 2(a)]. The flow field computed by numerical simulations by Liu et al.³⁰ revealed that the velocity profile first decreases when the fluid enters each element and subsequently increases, leading to the break of particle aggregates. A quantification of the breaking efficiency was, however, not reported. For the 20 μm particles, we observed several aggregates at the device inlet, i.e., before entering the series of trapezoidal elements [left panel of Fig. 5(a)]; at the end of the series, we observed a decrease of particle aggregates [middle and right panels of Fig. 5(a)]. We quantified the breaking efficiency of the designed microfluidic structure for different imposed pressure drops $ΔP$ and two different particle bulk concentrations, namely, $ϕb=0.2$ wt. % and $ϕb=0.3$ wt. %. For $ϕb=0.2$ wt. %, at the device inlet [Fig. 5(b)], isolated particles accounted for 60%–80% of the overall number, depending on the imposed pressure drop. An increase of the pressure drop was found to increase the number of aggregates, in line with the fact that a larger number of particles would pass from a large tube to a contraction, thus increasing the chances of aggregate formation.³⁸ While moving toward the outlet of the microfluidic device, we observed a drastic reduction of the number of aggregates [Figs. 5(c) and 5(d)], with more than 90% of the particles being isolated. For $ϕb=0.3$ wt. % [Figs. 5(e)–5(g)], we observed a similar trend, with the only difference that the larger particle concentration resulted in an expected increased number of aggregates at the inlet. Regardless of the concentration, around 90% of particles were isolated at the end of the device. We also performed the same study on particles with a diameter of 45 μm [Figs. 5(h)–5(k)], finding a similar trend with an even higher percentage of isolated particles (up to 95%).

In summary, we experimentally demonstrated that the array of trapezoidal elements can aid the breaking of particle aggregates regardless of the particle size and the imposed pressure drop. We observed that around 90%–95% of the overall particles were isolated at the end of the microfluidic device.

We now proceed to present the results on the particle ordering at different distances from the inlet of a straight capillary. Since particle aggregates have a detrimental effect on particle ordering,³² we employed the device described before followed by a straight glass capillary [Fig. 2(b)]. In this way, we reduced the number of particle aggregates to around 5% of the overall number, thus being able to evaluate the particle self-assembly dynamics with a reduced disturbance caused by the aggregates. For 20 μm particles, corresponding to a confinement ratio $β=0.2$⁠, particles at 5 cm from the capillary inlet (right after the end of the array of trapezoidal elements) were focused on the centerline but were not ordered [black rectangle in Fig. 6(a)], while, toward the end of the capillary, we observed a preferential spacing between consecutive particles [magenta rectangle in Fig. 6(a)], in agreement with previous observations.^13,31 We also quantified the dynamics of particle ordering at different distances from the capillary inlet and pressure drops. For an imposed pressure $ΔP=300$ mbar (⁠$Q≈$ 0.5 μL/min, $De≈13$⁠), particles at $L/D=500$ (5 cm from the channel inlet) did not display any ordering, as clear from the lack of a peak in the $S*$ distribution [Fig. 6(b)]. When increasing the distance from the inlet up to $L/D=2500$ (25 cm from the capillary inlet), we observed a clear peak for $S*≃5$⁠, in agreement with previous experimental³¹ and numerical³² results. The same trend was observed when increasing the pressure drop up to 600 mbar (⁠$Q≈$ 5.5 μL/min, $De≈146$⁠) [Figs. 6(c)–6(e)], in very good agreement with the results from the numerical simulations reported in Fig. 4. Similar results were also observed when the particle bulk concentration increased to $ϕb=0.3$ wt. % [Figs. 6(f)–6(i)].

We then repeated the same set of experiments for the 45 μm particles, corresponding to $β=0.45$ [Figs. 6(j)–6(n)]. At variance with the 20 μm particles, we observed a preferential spacing already at $L/D=1000$ (10 cm from the capillary inlet), confirming that a larger confinement ratio promoted a faster transition to a self-ordered structure, in good agreement with the numerical simulation results previously discussed. We also observed that the spatial distribution of $S*$ varied less significantly with the channel length compared to $β=0.2$⁠. For instance, at $ΔP=900$ mbar and $ΔP=1000$ mbar, there is no significant difference of the distribution of the particle spacing between $L/D=2000$ and $L/D=2500$⁠, compatible with the fact that the ordering dynamics tends to slow down after achieving a Gaussian-like shape around the $S*$ peak. This last experimental observation is again in agreement with the numerical simulation results displayed in Figs. 4 and S3–S6.

The experimental results just described represent a qualitative validation of the numerical simulations presented at the beginning of the previous section. Even though the constitutive equation employed for the numerical simulations does not quantitatively describe the rheological properties of the XG solution, both the model and the experimental fluid presented elastic and shear-thinning properties, thus making us highlight a few points that characterize the self-assembly of particle trains in viscoelastic shear-thinning solutions. Particle trains are observed at the centerline, in agreement with the fact that particles suspended in shear-thinning liquids are aligned at the channel centerline.^39,40 The transition from simple focusing to a self-ordered structure requires a certain channel length depending on the confinement ratio. This has been observed both in numerical simulations and in experiments, where the self-ordered structures formed closer to the channel inlet at $β≥0.4$ as compared to $β=0.2$⁠. By observing the histograms in Figs. 4 and 6, it is clear that the volumetric flow rate does not have a significant effect on the preferential spacing attained by the particles, which is controlled by the geometry of the channel and the particle linear concentration, in agreement with previous works.^29,31 An increase in the linear concentration led to an increase of particle doublets (Figs. 4, S5, and S6), also in agreement with recent experiments.³¹ Particle aggregation can be explained in light of the numerical simulations of Fig. 3, where three aligned particles experienced attractive forces below a critical distance. In other words, when the linear concentration increases, the number of particles approaching the channel grows, thus leading to the possibility of having particles closer to each other, where attractive forces are favored over repulsive ones.^29,31,32 The distance from the channel inlet remains another important parameter that affects the development of a ordered structure. The present numerical simulations (Fig. 4) and experiments (Fig. 6) clearly showed that viscoelastic ordering is a slow phenomenon and long distances from the inlet are required to attain a self-ordered structure. Furthermore, once the distribution of the interparticle distances approaches a Gaussian-like shape (with the exclusion of the particle aggregates), the ordering dynamics slows down significantly, with the distributions remaining substantially unaltered. In other words, a train of equally spaced particles in viscoelastic liquids may require a significant channel length, even larger than $L/D=3000$⁠, and such a distance increases for small confinement ratios. In this regard, it should also be mentioned that the migration dynamics of the particles toward the centerline in shear-thinning liquids that takes place before ordering is itself affected by the confinement ratio.⁴¹ Specifically, the migration velocity in a pressure-driven channel flow scales with $β2$⁠, requiring a longer channel also in the focusing region as the confinement ratio decreases.

Despite the quantitative effect of the confinement ratio on the required channel length to achieve ordering, the distributions of the interparticle distance are not significantly affected for linear concentrations up to $ϕ=0.3$⁠. As visible from the histograms in Figs. S3 and S4 (corresponding to $ϕ=0.1$ and $ϕ=0.2$⁠), self-similar distributions are found at decreasing distances from the channel inlet as the confinement ratio increases. As previously remarked, this behavior is a direct consequence of the scaling of the relative velocities with the confinement ratio that holds for sufficiently large interparticle distances. As the system is more crowded and the particles get closer, the scaling is lost and the magnitude of the scaled relative velocities is larger at small confinement ratios [see the left part of the panels in Fig. 3(b)]. Indeed, as visible in Fig. S3 corresponding to $ϕ=0.3$⁠, the distributions at a low confinement ratio show a higher number of particles in contact as compared to more confined geometries at the same “equivalent distance.”

It is worth mentioning that the self-similar dynamics with varying the confinement ratio is not observed in inertial ordering. Recent numerical simulations on staggered particle trains have reported a strong influence of the confinement ratio on the dynamics and stability of the train.²⁴ Specifically, unstable trains are observed for low Reynolds numbers and high confinement ratios due to the dominant interaction between particles over inertial effects. Similar instabilities are also found for small confinement ratios and high Reynolds numbers where the particle spacing fluctuates irregularly with large amplitude oscillations. Between these two limiting behaviors, for a value of the Reynolds number lower than a threshold, the staggered train is stable; otherwise, the interparticle spacing shows periodic oscillations without reaching a steady-state condition. The threshold increases with increasing the confinement ratio, i.e., higher confinement ratios tend to stabilize the train.

The formation of linear particle trains through inertial effects is also affected by the confinement ratio as recently reported by numerical simulations.²⁵ The spacing between the particles increases as the confinement ratio increases. This is a direct consequence of the displacement of the lateral equilibrium position attained by a single particle as the confinement ratio is varied. Specifically, as $β≪1$⁠, the particles approach the Segré and Silberberg equilibrium position, whereas, as the confinement ratio increases, the equilibrium position shifts toward the centerline due to steric effects,^42–45 affecting, in turn, the particle–wall hydrodynamic interactions and the average train velocity. Consequently, the required channel length to achieve ordering and the resulting particle spacing are a non-trivial function of the confinement ratio. In contrast, as fluid elasticity and shear-thinning “trap”' the particles at the channel centerline regardless of their size^39,40 for $β≥0.2$⁠, the ordering dynamics is roughly the same as the confinement ratio is varied, allowing to easily estimate the required channel length to achieve a certain ordering efficiency.

In this work, we performed numerical simulations and experiments to clarify the effect of the confinement ratio on the self-assembly dynamics of particle trains in viscoelastic liquids flowing in a straight microchannel. Both simulations and experiments provided a unified picture of the self-ordering dynamics that can be summarized according to the following results. For a fixed tube diameter, larger particles tended to self-order on shorter channel distances compared to smaller particles. The ordering dynamics was strongly controlled by the distance between consecutive particles, with a critical distance marking the transition between attractive and repulsive forces. Large particle concentrations resulted in the formation of particle aggregates that may suppress the train formation; to address this problem, we designed a microfluidic device made of an array of trapezoidal elements similar to Liu et al.³⁰ and demonstrated that such structure was able to significantly reduce the percentage of particle aggregates. The relative velocities between consecutive particles can be scaled with the fourth power of the confinement ratio when plotted as a function of the distance between the particle surfaces normalized by the channel diameter. The scaling holds for an interparticle spacing higher than a threshold, i.e., at relatively low linear concentrations (up to 0.2). Hence, in this case, the train dynamics at different confinement ratios is self-similar, leading to an evolution of the interparticle spacing that is simply rescaled with respect to the distance from the channel inlet. The experimental observations are consistent with the numerical predictions, and our findings are expected to hold for other polymer solutions presenting both elastic and shear-thinning properties.

Future works should look at clarifying the impact of attractive and repulsive forces between consecutive particles as a function of the fluid rheology, as only limited experiments exist on the subject.²⁹ Additional work is required to understand the potential effect of channel geometry (e.g., different cross sections and overall design) on the self-ordering dynamics. Finally, another important fluid dynamics problem is the stability of the train when subjected to flow disturbance. To the best of our knowledge, only a simple argument regarding train stability was previously introduced:²⁹ numerical simulations and experiments may look at how a fully formed particle train can react when subjected to flow disturbances such as obstacles or changes in the channel cross section.

See the supplementary material for (i) additional figures related to the relative velocity values between consecutive particles, similar to Fig. 3 and (ii) additional histograms of interparticle distance distribution, similar to Fig. 4.

F.D.G. acknowledges support from EPSRC Grant No. EP/S036490/1.

The authors have no conflicts to disclose.

F.D.G. and G.D.A. designed the study. A.J. performed the experiment and G.D.A. performed the numerical simulations. All the authors analyzed the data, wrote the manuscript, and approved the final submission.

The data that support the findings of this study are available from the corresponding author upon reasonable request.