4.1. Pairwise interactions (PIs)

The simulation method for PIs involves calculating the particle trajectories of two SAPs initially separated widely in the flow-direction, $1\ll |\Delta r_{0,1}|$ . The relative initial positions in the gradient, $\Delta r_{0,2}=O(1)$ , and vorticity, $\Delta r_{0,3}=O(1)$ , directions are randomly chosen in a specified region of area $\mathbb{S} =O(1)$ . As the two SAPs approach each other due to the relative motion induced by the background SSF, both their translational and rotational dynamics are impacted by long-range HIs and short-range mechanical contacts. These interactions modify the particle stresslet relative to its isolated-particle value, described in § 2, during the course of the trajectory, even though the particles ultimately return to their stable equilibrium orientation $\boldsymbol{p}_S$ at the end of the interaction.

These transient stresses induced during the interaction provide $O(n^2)$ corrections to the particle stresses in the suspension. These transient stresses for SAPs, can be separated into three distinct categories: (i) $ \boldsymbol{S}_p$ , the stresslet arising from the transient deviation of the particle orientation from $\boldsymbol{p}_S$ ; (ii) $\boldsymbol{S}_f$ , the stresslet associated with the velocity disturbance from long-range HIs incorporated through $\boldsymbol{f}_{\textit{HI}}$ ; and (iii) $ \boldsymbol{S}_{\textit{con}}$ , the stresslet resulting from contact forces that prevent particles from passing through each other. The particle stress $\boldsymbol{\sigma }_p$ in the suspension can be described as

(4.3) \begin{align} \boldsymbol{\sigma }_p = n\langle \boldsymbol{S}_{\textit{Hyd}}\rangle _J + n^2 \left (\frac {\langle \boldsymbol{S}_p\rangle }{n} + \frac {\langle \boldsymbol{S}_f\rangle }{n} + \frac {\langle \boldsymbol{S}_{\textit{con}} \rangle }{n}\right ), \end{align}

where the $O(n^2)$ terms are the ensemble average values of the three different contributions described earlier, $\langle \boldsymbol{\cdot }\rangle$ representing the ensemble average over all realizations of PIs. This decomposition enables a clear separation of orientational, hydrodynamic and contact-driven contributions to the rheology.

The value of $\langle \boldsymbol{S}_p\rangle$ can be evaluated by computing the change in $\boldsymbol{S}_{\textit{Hyd}}$ relative to $\langle \boldsymbol{S}_{\textit{Hyd}}\rangle _J$ , during the trajectory. The ensemble average of $\boldsymbol{S}_{p}=(\boldsymbol{S}_{\textit{Hyd}} - \langle \boldsymbol{S}_{\textit{Hyd}}\rangle _J)$ is

(4.4)

where the inner time integral is taken over each PI trajectory, and the outer integral is taken over all values of the initial cross-streamline separations of the two particles in the gradient–vorticity plane. We used a Monte Carlo summation to evaluate the integral because it is computationally more efficient for tumbling particles, which have four additional degrees of freedom corresponding to initial orientations of the two particles (Borker & Koch Reference Borker and Koch2023). This orientation distribution is described by a Jeffery orbit distribution and a phase angle chosen uniformly between $(-\pi ,\pi ]$ . The details of obtaining the orbit distribution are described in section II.D of Borker & Koch (Reference Borker and Koch2023). The value of $\langle \boldsymbol{S}_p\rangle$ can be alternatively computed as

(4.5) \begin{align} \langle \boldsymbol{S}_p\rangle = n \frac {\mathbb{S}}{N} \sum _{j=1}^{N} \Delta r_{0,2,j} \int ^{\infty }_{-\infty } {\rm d}t \left (\boldsymbol{S}_{\textit{Hyd}} - \langle \boldsymbol{S}_{\textit{Hyd}}\rangle _J \right ), \end{align}

where $N$ is the number of different initial relative positions of two SAPs in the gradient–vorticity plane which are chosen uniformly over a specified region of area $\mathbb{S}=O(1)$ , and $\Delta r_{0,2,j}$ is the initial relative separation in the gradient direction that weights the contributions from different trajectories based on the SSF driven particle flux across the gradient–vorticity plane (Da Cunha & Hinch Reference Da Cunha and Hinch1996; Borker & Koch Reference Borker and Koch2023).

Direct evaluation of the integrals in (4.4)–(4.5) leads to non-convergent integrals because the integrand in (4.5), $\boldsymbol{S}_{\textit{Hyd}} - \langle \boldsymbol{S}_{\textit{Hyd}}\rangle _J$ , is proportional to the velocity disturbance generated by a neighbouring ring which decays as $O(\Delta r^{-3})$ , $\Delta r$ being the interparticle distance. To obtain a physically meaningful, absolutely convergent integral, we therefore need to adopt a renormalization strategy that subtracts a contribution whose ensemble average is identically zero but whose far-field functional form matches that of the original integrand at $O(\Delta r^{-3})$ . This approach follows the classical renormalization framework introduced by Batchelor & Green (Reference Batchelor and Green1972) for spheres and subsequently extended to non-spherical particles by Borker & Koch (Reference Borker and Koch2023). The motivation for the term used for the renormalization procedure can be determined by understanding the fact that PIs occur in a dilute suspension of perfectly aligned SAPs rather than in a background Newtonian fluid. Consider, the linearized change in a SAP’s orientation, $\boldsymbol{\Delta p}$ , due to the velocity disturbance generated by a uniform suspension of perfectly aligned SAPs. The ensemble average value of $\boldsymbol{\Delta p}$ should be zero, because the mean velocity gradient in a uniform suspension of perfectly aligned SAPs is the same as that in the imposed shear flow, i.e. $\boldsymbol{\varOmega }_\infty +\boldsymbol{E}_\infty$ . This ensemble average can also be evaluated using pairwise trajectory approach by computing $\boldsymbol{\Delta p}$ for a SAP as it is influenced by all possible approaches of perfectly aligned SAPs.

The term to be subtracted from (4.5) for the renormalization procedure is

(4.6) \begin{align} \boldsymbol{\Delta S}^{\textit{RN}}_{\textit{Hyd}}(\boldsymbol{\Delta p},\boldsymbol{p}_S) & = \xi _2 \bigg (\boldsymbol{E}_\infty \boldsymbol{\cdot }\left (\boldsymbol{\Delta p}\,\boldsymbol{p}_S + \boldsymbol{p}_S\,\boldsymbol{\Delta p} \right )+\left (\boldsymbol{\Delta p}\,\boldsymbol{p}_S + \boldsymbol{p}_S\,\boldsymbol{\Delta p} \right )\boldsymbol{\cdot }\boldsymbol{E}_\infty \nonumber \\& \quad - \frac {2}{3}\boldsymbol{E}_\infty :\left (\boldsymbol{\Delta p}\,\boldsymbol{p}_S + \boldsymbol{p}_S\,\boldsymbol{\Delta p}\right )\boldsymbol{I}\bigg ) \nonumber \\& \quad + \xi _3 \bigg (\boldsymbol{E}_\infty : \big (\boldsymbol{\Delta p}\,\boldsymbol{p}_S \boldsymbol{p}_S \boldsymbol{p}_S + \boldsymbol{p}_S\,\boldsymbol{\Delta p} \, \boldsymbol{p}_S \boldsymbol{p}_S +\boldsymbol{p}_S \boldsymbol{p}_S \, \boldsymbol{\Delta p}\, \boldsymbol{p}_S\nonumber \\& \quad + \boldsymbol{p}_S \boldsymbol{p}_S \boldsymbol{p}_S \,\boldsymbol{\Delta p}\big ) - \frac {1}{3}\boldsymbol{E}_\infty :\left (\boldsymbol{\Delta p}\,\boldsymbol{p}_S+\boldsymbol{p}_S\,\boldsymbol{\Delta p}\right )+ 2\boldsymbol{p}_S \boldsymbol{p}_S \, \boldsymbol{\Delta p}\cdot\boldsymbol{p}_S\boldsymbol{I}\bigg ), \end{align}

where the rate of change of $\Delta \boldsymbol{p}$ with time is given by

(4.7) \begin{eqnarray} \frac {{\rm d}\boldsymbol{\Delta p}}{{\rm d}t} & = & \boldsymbol{\omega }_{RN}\times \boldsymbol{p}_S + \boldsymbol{\Delta p}\boldsymbol{\cdot }\boldsymbol{\varOmega }_\infty + \lambda \big (\boldsymbol{E}_\infty \boldsymbol{\cdot }\boldsymbol{\Delta p} \nonumber \\ & - & \boldsymbol{E}_\infty :\left (\boldsymbol{\Delta p}\,\boldsymbol{p}_S \boldsymbol{p}_S + \boldsymbol{p}_S \boldsymbol{\Delta p}\,\boldsymbol{p}_S + \boldsymbol{p}_S \boldsymbol{p}_S \, \boldsymbol{\Delta p}\right )\big ). \end{eqnarray}

Here, $\boldsymbol{\omega }_{RN}$ is the angular velocity induced from the velocity disturbance of a perfectly aligned SAP translating along the fluid streamline with the same initial cross-streamline separation as the original trajectory calculation. Using the same initial cross-streamline separations as the original calculation reduced the statistical error associated with the Monte Carlo summation (Borker & Koch Reference Borker and Koch2023). The first term is rotation caused by the gradient in the velocity disturbance produced by a suspension of non-interacting SAPs in the form of $\boldsymbol{\omega }_{RN}$ . The second and third terms together are the Jeffery rotation rate correct to linear order in $\boldsymbol{\Delta p}$ (obtained through a Taylor series expansion of (2.1)) and is induced by the deviation of the particle’s orientation from $\boldsymbol{p}_S$ . For a given set of initial conditions, the value of ${\rm d}\boldsymbol{\Delta p}/{\rm d}t$ for the renormalized equation matches the original calculation at $O(\Delta r^{-3})$ (Borker & Koch Reference Borker and Koch2023). This implies that $\boldsymbol{\Delta p}+ \boldsymbol{p}_S$ from the linearized equation also matches $\boldsymbol{p}$ from the original calculation for each trajectory at $O(\Delta r^{-3})$ . Therefore, $\boldsymbol{S}_{\textit{Hyd}} - \langle \boldsymbol{S}_{\textit{Hyd}}\rangle _J$ also matches to $O(\Delta r^{-3})$ . The renormalized integral for computing $\langle \boldsymbol{S}_p\rangle$ is

(4.8) \begin{align} \langle \boldsymbol{S}_p\rangle = n \frac {\mathbb{S}}{N} \sum _{j=1}^{N} \Delta r_{0,2,j} \int ^{\infty }_{-\infty } {\rm d}t \left (\boldsymbol{S}_{\textit{Hyd}} - \langle \boldsymbol{S}_{\textit{Hyd}}\rangle _J - \boldsymbol{\Delta S}^{\textit{RN}}_{\textit{Hyd}} \big(\boldsymbol{\Delta p},\boldsymbol{p}_S \big)\right ). \end{align}

Figure 8(a) illustrates the original (dashed black line) and renormalized (solid red line) orientation trajectory of the particles for a given initial cross-streamline separation of the two particles which is much larger than unity. In the renormalized trajectory, the fluctuations in the particle orientation are nearly zero within the accuracy of the numerical simulation. Figure 8(b) shows the contribution of the $12$ component of $\boldsymbol{S}_p$ from original and renormalized ( ${S}_{p,12}-\Delta S^{\textit{RN}}_{Hyd,12}$ ) calculations, which also depicts a similar trend. This suggests that the integral in (4.8) rapidly converges with increasing values of the initial cross-streamline separation of the particles. For tumbling particles, the constant orientation vector $\boldsymbol{p}_S$ in (4.6)–(4.8), is replaced by a time varying Jeffery trajectory of the particle described by (1.2), and the details are available in section II.C of Borker & Koch (Reference Borker and Koch2023).

Figure 8.

Renormalization of PI trajectories for SAPs at large cross-streamline separations. (a) The time evolution of the flow component of the orientation vector of a SAP relative to its equilibrium value for the original calculation, $({p}_1-{p}_{S,1})/\delta$ (dashed line) and the renormalized calculation, $({p}_1-{p}_{S,1} - {\Delta p}_1 )/\delta$ . Time evolution of the 12 components of the stresslet tensor, (b) $\boldsymbol{S}_{p}$ and (c) $\boldsymbol{S}_f$ . The cross-streamline separation of the particles was $(\Delta r_{0,2},\Delta r_{0,3}) = (6.16, 1.09)$ .

The HI contribution to the stresslet $\langle \boldsymbol{S}_f\rangle$ is given by

(4.9) \begin{align} \langle \boldsymbol{S}_f\rangle = n \frac { \mathbb{S}}{N} \sum _{j=1}^{N} \Delta r_{0,2,j} \int ^{\infty }_{-\infty } {\rm d}t \left ( \left [\oint ds\left (\frac {\boldsymbol{r}_{\!c} \boldsymbol{f}_{\textit{HI}}+\boldsymbol{f}_{\textit{HI}}\boldsymbol{r}_{\!c}}{2}-\frac {\boldsymbol{f}_{\textit{HI}}\boldsymbol{\cdot }\boldsymbol{r}_{\!c}}{3}\right ) \right ]- \boldsymbol{S}^{\textit{RN}}_{\kern-1pt f} \right ), \end{align}

where $\boldsymbol{r}_{\!c}$ is the centreline position of the ring cross-section, $s$ denotes the position along $\boldsymbol{r}_{\!c}$ , $ds$ is the elemental length along $\boldsymbol{r}_{\!c}$ and the inner integral is evaluated along $\boldsymbol{r}_{\!c}$ . The renormalization HI stresslet $\boldsymbol{S}^{\textit{RN}}_{\kern-1pt f}$ is the stresslet associated with the force per unit length $\boldsymbol{f}_{\textit{HI}}^{\textit{RN}}$ on a SAP due to the velocity disturbance produced by a non-interacting SAP aligned at their equilibrium orientation. For each initial cross-streamline separation, we evaluate $\boldsymbol{f}_{\textit{HI}}^{\textit{RN}}$ due to the velocity disturbance from a SAP translating with the SSF and has a fixed orientation corresponding to $\boldsymbol{p}=\boldsymbol{p}_S$ (Borker & Koch Reference Borker and Koch2023). The contributions from $\boldsymbol{S}_f$ and $\boldsymbol{S}^{\textit{RN}}_{\kern-1pt f}$ are identical at $O(\Delta r^{-3})$ when $\Delta r\gg 1$ , where $\boldsymbol{S}_f$ is the term inside the brackets $[ \boldsymbol{\cdot }]$ in (4.9). For a trajectory with a large cross-streamline separation, the time evolution of ${S}_{f,12}$ and ${S}_{f,12}-{S}^{\textit{RN}}_{f,12}$ are depicted in figure 8(c), indicating that the renormalized integral converges rapidly. For tumbling rings, the implementation is described in section II.C of Borker & Koch (Reference Borker and Koch2023). The ensemble average value of $\boldsymbol{S}^{\textit{RN}}_{\kern-1pt f}$ is identically zero as the ensemble average of the velocity disturbance across the suspension is identically zero. This renormalization method of obtaining an absolutely convergent summation is similar to the one used for spheres (Batchelor & Green Reference Batchelor and Green1972).

The contribution of the normal contact forces to the stresslet is given by

(4.10) \begin{align} \langle \boldsymbol{S}_{\textit{con}}\rangle = n\frac { \mathbb{S}}{N} \sum _{j=1}^{N} \Delta r_{0,2,j} \int ^{\infty }_{-\infty } {\rm d}t \frac {1}{2}\left (\frac {\boldsymbol{\Delta r} \boldsymbol{F}_{\textit{con}}+\boldsymbol{F}_{\textit{con}}\boldsymbol{\Delta r}}{2}-\frac {\boldsymbol{F}_{\textit{con}}\boldsymbol{\cdot }\boldsymbol{\Delta r}\boldsymbol{I}}{3}\right ), \end{align}

where $\boldsymbol{F}_{\textit{con}}$ is the contact force needed to prevent particles from passing through each other during the interaction (see section II.A of Borker & Koch (Reference Borker and Koch2023)).

Figure 9 depicts the transient trajectories of pertinent variables of interest of an $A=40$ SAP during three distinct PIs. These trajectories have some of the largest contributions to the $O(n^2)$ particle stresses induced by PIs. Figure 9(a) illustrates the translational trajectories of the SAP during these PIs. The trajectory is asymmetric about the gradient–vorticity plane and the relative cross-streamline separation of the particles also changes during the interaction. This cross-streamline displacement is responsible for inducing translational hydrodynamic diffusion of particles in the suspension (Borker & Koch Reference Borker and Koch2023). This translational hydrodynamic diffusivity can be quantified from the ensemble average of the cross-streamline displacements (Da Cunha & Hinch Reference Da Cunha and Hinch1996; Borker & Koch Reference Borker and Koch2023)

(4.11) \begin{align} D_{ii} = \frac {n}{2}\frac {\mathbb{S}}{N}\sum _{j=1}^{N} \Delta r_{0,2,j} \left (\Delta r_{f,i}\right )^2, \end{align}

where $\Delta r_{f,i}$ , for $i\in \{2,3\}$ , is the cross-streamline displacement of the particle from its initial streamline. Figure 19 depicts the hydrodynamic diffusivities, $D_{22}$ and $D_{33}$ , of T-rings and tori with varying aspect ratios, indicating that SAPs have orders of magnitude smaller dispersion in a suspension compared with tumbling particles. This is primarily caused by the $O(1)$ cross-streamline displacements caused by mechanical contacts between tumbling particles (Borker & Koch Reference Borker and Koch2023) causing $O(1)$ cross-streamline displacements, while SAPs mainly interact through weak HIs (Borker et al. Reference Borker, Stroock and Koch2024).

Figure 9(b) shows the $O(\delta )$ displacement of the particle’s orientation relative to $\boldsymbol{p}_S$ that illustrates the wobbling of a SAP around its stable orientation during PIs with another SAP. The magnitude of the wobbling was found to be $O(\delta )$ for trajectories which had large contributions to the suspension properties. Figures 9(c) and 9(d) represent the corresponding changes in $S_{p,12}$ and $S_{f,12}$ as a function of the separation between the centres of mass of the two SAPs in the flow direction. These results demonstrate that even when particles approach each other, the induced stresses remain small in the range $10^{-3}$ – $10^{-2}$ .

Figure 9.

Trajectories of SAPs due to PIs with another SAP that have some of the largest contributions to the particle stresses. (a) Relative separation of the two SAPs in the flow and gradient directions. (b) Change in the particle’s orientation vector projected in the flow–vorticity plane. The change in (c) $S_{p,12}$ and (d) $S_{f,12}$ as a function of the separation between the centre of masses of the two SAPs in the flow direction. The symbol o in the graphs represents the point of minimum separation of the particles’s centre of mass in the flow direction. Different lines indicate three different initial cross-streamline separations.

The components of the stresslet tensor that contribute to the shear viscosity of the suspension are denoted by $\langle S_i\rangle =\langle \boldsymbol{S}_i\rangle :\hat {\boldsymbol{e}}_1\hat {\boldsymbol{e}}_2$ , with $i\in \{p,f,con\}$ . Figure 10 shows the variation of $\langle {S}_i\rangle$ with the particle aspect ratio $A$ . It is evident that shear stresses associated with PIs are substantially smaller for SAPs compared with ETs. This difference arises primarily from the contrasting nature of the interactions. For ETs, the dominant contributions come from collisions between a tumbling and a temporarily aligned particles (Rahnama et al. Reference Rahnama, Koch and Shaqfeh1995; Borker & Koch Reference Borker and Koch2023), which cause large $O(1)$ deviations in their positions and orientations, as illustrated in figure 4(a–d) of Borker & Koch (Reference Borker and Koch2023). Hydrodynamic interactions would cause smaller $O(\epsilon )$ changes in the particle’s orientation, position and stresslet. Consequently, collisions have the dominant contribution to particle stresses associated with PIs of ETs (Borker & Koch Reference Borker and Koch2023). On the contrary, SAPs mainly interact through purely HIs rather than mechanical contacts. This is reflected in their much smaller magnitude of contact stresses compared with ETs, as depicted in figure 10(c), and the lower collision frequency as shown in figure 8(a) of Borker et al. (Reference Borker, Stroock and Koch2024). The trajectories in figure 9 also illustrate the weak interaction of SAPs. The permanent alignment of SAPs near the fluid lamellae not only reduces the leading-order shear stresses compared with ETs but also lowers the strength of interparticle interactions and the stresses generated during these interactions.

Figure 10.

Contributions to the shear stress from PI. Variation of (a) $\langle {S}_p\rangle$ , the contribution driven by changes in orientation from the isolated particle trajectory (for tumblers) or the equilibrium orientation (for SAPs); (b) $\langle {S}_f\rangle$ , the contribution from pairwise HIs; and (c) $\langle {S}_{\textit{con}}\rangle$ , the contribution from mechanical contacts; with the aspect ratio, $A$ . Lines represent SBT scaling with aspect ratio summarized in table 3.

We provide further insights from scaling arguments into the origin of the smaller stresses from interparticle interactions in SAPs compared with ETs. Our simulations indicated that the dominant PIs occurred at $O(1)$ gradient separations, implying that the duration of interaction is also expected to be $O(1)$ . The scaling of $\langle {S}_i\rangle$ would thus be equivalent to the scaling of the integrand in (4.8)–(4.10). The variation of $\langle {S}_p\rangle$ with $A$ is depicted in figure 10(a) for both SAPs and ETs, and the behaviour mirrors the trends observed in the aspect ratio dependence of $\langle S_{\textit{Hyd}}\rangle _J$ illustrated in figure 3(a). The HIs of SAPs induce an $O(\epsilon \delta )=O(\epsilon \sqrt {\Delta A/A^*}/A)$ change in particle orientation relative to $\boldsymbol{p}_S$ , while contacting interactions only cause an $O(\delta ^2)=O(\Delta A/(A^*A^2))$ change. Therefore, $\langle {S}_p\rangle$ scales as $O(n\epsilon ^2\Delta A/(A^*A^2))$ . For ETs, by contrast, mechanical contact between a tumbling and a temporarily aligned ET, occurring with a frequency of $O(\delta _T)$ , induces a large $O(1)$ change in the orientation relative to the particle’s periodic Jeffery trajectory for each such interaction. This results in an overall contribution of $O(n\epsilon \delta _T)=O(n\epsilon /A)$ to $\langle {S}_p\rangle$ (Borker & Koch Reference Borker and Koch2023). The scaling also indicates that the magnitude of $\langle {S}_p\rangle$ for SAPs with aspect ratios that are sufficiently far away from $A^*$ , i.e. $A\gg A^*$ , is expected to be $O(1/A^*)\ll 1$ smaller than for ETs which is consistent with the observation in figure 10(a). The solid lines in figure 10(a) represents the aforementioned scalings for SAPs and ETs and their values are described in table 3.

The HI shear stresslet $\langle {S}_f\rangle$ is driven by the velocity disturbance from neighbouring particles, $|\boldsymbol{u}^\prime _\infty |$ , which is $O(\epsilon \delta )=O( \epsilon \sqrt {\Delta A/A^*}/A)$ , based on SBT (Borker & Koch 2019, Reference Borker and Koch2023). The SBT indicates the magnitude of the stresses generated from these velocity disturbances should be $O(\epsilon |\boldsymbol{u}^\prime _\infty |)$ and therefore $\langle {S}_f\rangle$ scales as $O(\epsilon ^2 \sqrt {\Delta A/A^*}/A)$ . In comparison, for tumbling particles, $\langle {S}_f\rangle =O(\epsilon ^2\delta _T)=O(\epsilon ^2/A)$ , exhibiting a similar aspect ratio scaling for $A\gg A^*$ , albeit with a numerically larger magnitude due to stronger PIs.

Finally, the contribution from contact forces for SAPs is also notably lower than for ETs, as depicted in figure 10(c). This difference arises because the $O(\delta ^2)=O(\Delta A/(A^*A^2))$ collision frequency of SAPs is much smaller compared with the $O(1/A)$ collision frequency of ETs, dominated by collisions between tumbling and temporarily aligned particles (Borker & Koch Reference Borker and Koch2023). Moreover, $\langle S_{\textit{con}}\rangle$ was found to be negative for SAPs. The negative sign can be explained by the observation during the typical binary collision events that the contact force was approximately oriented along the gradient direction $\pm \hat {\boldsymbol{e}}_2$ when the point of contact relative to the particle’s centre of mass was approximately oriented along $\mp \hat {\boldsymbol{e}}_1$ . This indicates that the contribution to the $\hat {\boldsymbol{e}}_1\hat {\boldsymbol{e}}_2$ component of $\langle \boldsymbol{S}_{\textit{con}}\rangle$ should have a negative sign. Our simulations also indicate that the magnitude of these contributions to the ensemble average were an order of magnitude smaller compared with $\langle S_{p}\rangle$ and $\langle S_{f}\rangle$ .

Overall, our dynamic simulations demonstrate the key physical insight that the intrinsic flow alignment of SAPs limits the strength of PIs between particles and any additional stresses generated in the process. We demonstrate that SAPs primarily interact through weak long range HIs unlike tumbling planar particles which interact through binary collisions. The SAPs with other geometries, described in Borker et al. (Reference Borker, Stroock and Koch2018, Reference Borker, Stroock and Koch2024) should also follow the same trends described here.

4.2. Tuneable microstructure and viscosity using SAP–ET mixtures

In addition to PIs between SAPs, we also studied the PIs between a SAP and an ET. Figure 11(a) describes the change in the position of an ET relative to the SAP in the flow–gradient plane. The trajectories in figure 11 represent cases that generate some of the largest stresses in a SAP due to interactions with an ET. The thicker sections of each trajectory represent portions of the trajectory when the particles were in mechanical contact with each other. Figure 11(b) shows the change in the particle’s orientation projected on the flow–gradient plane. This indicates that not all collisions with a tumbling particle cause a SAP to tumble. Notice that the orientational disturbances induced by the ET are significantly larger than those present during a SAP–SAP interaction. Figures 11(c) and 11(d) illustrate the change in $S_{p,12}$ and $S_{f,12}$ , respectively, as a function of the relative separation of the ET from the SAP in the flow direction. Firstly, the magnitudes of the stresses are much larger than those described in figures 9(c) and 9(d). The sharp rise is caused by both physical collisions (thicker segment of the trajectory) and long range HIs (thinner segments of the trajectory) induced by the passing tumbling particle. The tumbling particles induce minimal stresses when they are sufficiently far away from the SAP even within an $O(1)$ distance. The major contribution to the stress on the SAP was caused by its deviation from its equilibrium orientation, i.e. $S_{p,12}$ compared with $S_{f,12}$ or $S_{con,12}$ .

Figure 11.

Trajectories of $A=40$ SAPs due to PIs with $A=40$ torus that have some of the largest contributions to the particle stresses. (a) Separation of the ET relative to the SAP in the flow and gradient directions. (b) Change in the orientation vector of the SAP projected in the flow–gradient plane. The change in (c) $S_{p,12}$ and (d) $S_{f,12}$ as a function of the relative separation of the ET from the SAP in the flow direction.

Figure 12.

Properties for SAP–ET mixtures. (a) Specific viscosity and (b) flow-alignment parameter of a mixture of $A=40$ SAPs and ETs as a function of the non-dimensional number density, $n$ , at different particle fractions of ETs, $\chi$ . The shaded area corresponds to 95 % interval confidence intervals.

These simulations in conjunction with SAP–SAP interactions studied here and the ET–ET interactions of Borker & Koch (Reference Borker and Koch2023) allow one to estimate how the orientational microstructure and specific viscosity change with small additions of ETs in SAP suspensions. The specific viscosity of a suspension of SAPs can be precisely increased by adding a small and controlled fraction of ETs to the suspension. The specific viscosity, $\Delta \mu$ , for a mixture of SAPs and ETs can be computed from the following equation:

(4.12) \begin{eqnarray} \Delta \mu = &n&\left ((1-\chi )\langle S_{\textit{Hyd}}\rangle _{J,S} + \chi \langle S_{\textit{Hyd}}\rangle _{J,T}\right ) \nonumber \\ &+& n^2\left ((1-\chi )^2 \langle S\rangle _{\textit{PI},S}+ \chi (1-\chi )\langle S\rangle _{PI,ST} + \chi ^2\langle S\rangle _{\textit{PI},T} \right ), \end{eqnarray}

where $\chi \ll 1$ is the particle fraction of ETs, $\langle S\rangle _{PI}=\langle S_p\rangle +\langle S_f\rangle +\langle S_{\textit{con}}\rangle$ ; the subscripts $S,T$ correspond to values from PI trajectories of two SAPs and two ETs, respectively; and the subscript ST corresponds to values evaluated from PI trajectories between a SAP and an ET (see S4 of Borker et al. (Reference Borker, Stroock and Koch2024)). Equation (4.12) is valid for $\epsilon ^2\delta ^2/(\epsilon ^2 \delta ^2 + \delta _T)\ll \chi \ll 1$ , wherein the steady state distribution of orbit constants of the ETs is influenced mainly by interactions with other ETs rather than SAPs. Equation (4.12) suggests that $\Delta \mu$ is linearly proportional to $\chi$ at $O(n)$ and its value changes from the value corresponding to SAPs at $\chi =0$ to the value corresponding to ETs at $\chi =1$ . Adding a small controlled number of ETs to SAP suspensions provides a flexible route for independently tuning the suspension viscosity and orientational microstructure. Our simulations reveal that the dominant $O(n^2)$ contributions to the viscosity arise from cross-interactions between SAPs and ETs rather than from ET–ET interactions.

To illustrate the predictions of (4.12), figure 12(a) shows $\Delta \mu$ correct to $O(n^2)$ for a mixture of $A=40$ particles for three different particle fractions of ETs: $\chi =0.05,0.10, 0.2$ . Figure 12(a) indicates that the specific viscosity of a mixture with 5 % ETs, i.e. $\chi =0.05$ , is only approximately 14 % of the value of $\Delta \mu$ for a suspension of pure ETs at a non-dimensional number density of $n=0.1$ . The alignment of the particles relative to the flow–vorticity plane would also change with $n$ and $\chi$ , and this can be quantified using a mixture average of the flow-alignment parameter $k$ defined as $(3/2)(1-\langle p_2^2 \rangle )$ (Borker & Koch Reference Borker and Koch2023; Borker et al. Reference Borker, Stroock and Koch2024). The flow-alignment parameter $k$ is 0 when the centreline of the ring-shaped particles is perfectly aligned in the flow–vorticity plane and 1 when the particles have an isotropic orientation distribution. The variation of $k$ with $n$ for different mixtures of $A=40$ particles is shown in figure 12(b), where the values were obtained from the simulation results of Borker et al. (Reference Borker, Stroock and Koch2024) (see (S4.1) and section S4 in their work). Therefore, one could choose the appropriate value of the particle fraction of ETs, $\chi$ , at a given non-dimensional number density, $n$ , to obtain the desired microstructural anisotropy and shear viscosity using results described in figure 12. Furthermore, by dynamically controlling the fraction of ETs in the suspension one could also construct composites with precisely controlled gradients in its orientational microstructure, which could be of interest to material scientists.

4.3. Multiparticle interactions (MIs)

We also conducted dynamic MI simulations within a periodic box to compute the effect of MIs on the rheology of SAP suspensions at $n=O(1)$ . We followed the methodology established for rigid fibre suspensions (Butler & Shaqfeh Reference Butler and Shaqfeh2002; Tornberg & Shelley Reference Tornberg and Shelley2004; Saintillan, Darve & Shaqfeh Reference Saintillan, Darve and Shaqfeh2005; Snook et al. Reference Snook, Davidson, Butler, Pouliquen and Guazzelli2014; Butler & Snook Reference Butler and Snook2018). These simulations track changes in the stresslet, translational and rotational velocities of the particle induced by the $O(\epsilon \delta )$ velocity disturbances generated by the leading-order force per unit length of neighbouring particles and the contact forces over time. The methodology is explained in detail in Saintillan et al. (Reference Saintillan, Darve and Shaqfeh2005) and the exact implementation for rings is given in Borker et al. (Reference Borker, Stroock and Koch2024). The particle stresses are given by

(4.13) \begin{align} \boldsymbol{\sigma }_p = n\left (\langle \boldsymbol{S}_{\textit{Hyd}}\rangle _{\textit{MI}}+ \langle \boldsymbol{S}_f\rangle _{\textit{MI}}+ \langle \boldsymbol{S}_{\textit{con}} \rangle _{\textit{MI}}\right ), \end{align}

where $\langle \boldsymbol{S}_{\textit{Hyd}}\rangle _{\textit{MI}}$ corresponds to the stresses on each particle from the imposed SSF, $\langle \boldsymbol{S}_f\rangle _{\textit{MI}}$ is the stresslet associated with the velocity disturbances generated by the presence of particles in the suspension and $\langle \boldsymbol{S}_{\textit{con}} \rangle _{\textit{MI}}$ is the stresslet associated with mechanical contact events. The portion of the stresslet termed as $\langle \boldsymbol{S}_{\textit{Hyd}}\rangle _{\textit{MI}}$ is computed as the time average of the volume averaged value of $\boldsymbol{S}_{\textit{Hyd}}$ in the suspension and given by

(4.14) \begin{eqnarray} \langle \boldsymbol{S}_{\textit{Hyd}}\rangle _{\textit{MI}} = \frac {1}{T_{Me}} \int _{T_e}^{T_M} {\rm d}t \langle \boldsymbol{S}_{\textit{Hyd}}\rangle _t =\frac {1}{T_{Me}} \int _{T_e}^{T_M} {\rm d}t \frac {1}{N_p} \sum _{j=1}^{N_p} \boldsymbol{S}_{\textit{Hyd}}(\boldsymbol{p}_j), \end{eqnarray}

where $T_M$ is the time over which the simulation is carried out, $T_e$ is the equilibration time for the particles in the periodic box to reach their statistical steady state orientation distribution, $ T_{Me} = T_M - T_e$ , $N_p$ is the number of particles in the periodic box, $j$ represents the index of the particle in the periodic box and $\langle \boldsymbol{S}_{\textit{Hyd}}\rangle _t$ is the transient volume average value of $\boldsymbol{S}_{\textit{Hyd}}$ in the suspension. Our simulations indicated that the suspension remained homogenous during the course of the simulation. This enables us to represent the volume averaged particle properties as an average over all the particles in the periodic box as described in (4.14). The time averaging was carried out for durations that were much greater than $O(\delta ^{-1})=O(A\sqrt {A^*/\Delta A})$ to ensure that the orientational distribution of SAPs relaxed to their statistical steady state distribution (see S2 of the Supplementary material of Borker et al. (Reference Borker, Stroock and Koch2024) for additional details). The hydrodynamic stresslet $\langle \boldsymbol{S}_{\textit{Hyd}}\rangle _{\textit{MI}}$ is equivalent to $\langle \boldsymbol{S}_{\textit{Hyd}}\rangle _{J}+\langle \boldsymbol{S}_{p}\rangle$ from PIs when $n\ll 1$ .

Similarly, $\langle \boldsymbol{S}_{f}\rangle _{\textit{MI}}$ is calculated as the time average of the volume averaged value of the stresslet associated with the velocity disturbances generated by all neighbouring particles and given by

(4.15) \begin{eqnarray} \langle \boldsymbol{S}_f\rangle _{\textit{MI}} &=& \frac {1}{T_{Me}} \int _{T_e}^{T_M} {\rm d}t \langle \boldsymbol{S}_{f}\rangle _t \nonumber \\ &=&\frac {1}{T_{Me}} \int _{T_e}^{T_M} {\rm d}t\frac {1}{N_p} \sum _{j=1}^{N_p} \oint _{\mathbb{S}_j} ds\left (\frac {\boldsymbol{r}^j_{c} \boldsymbol{f}^j_{\textit{HI}}+\boldsymbol{f}^j_{\textit{HI}}\boldsymbol{r}^j_{c}}{2}-\frac {\boldsymbol{f}^j_{\textit{HI}}\boldsymbol{\cdot }\boldsymbol{r}^j_{c}}{3}\right ), \end{eqnarray}

where the line integral is carried out over the centreline of the jth particle, $\mathbb{S}_j$ ; $\boldsymbol{f}^j_{\textit{HI}}$ is the force per unit length on the centreline of the jth particle associated with the velocity disturbances from all other particles in the suspension; $\boldsymbol{r}^j_{c}$ is the centreline position of the jth ring; and $\langle \boldsymbol{S}_{f}\rangle _t$ is the transient volume average value of the stresslet associated with velocity disturbances in the suspension. The stresslet from contact stresses, $\langle \boldsymbol{S}_{\textit{con}}\rangle _{\textit{MI}}$ , is given by

(4.16) \begin{align} \langle \boldsymbol{S}_{\textit{con}}\rangle _{\textit{MI}} & = \frac {1}{T_{Me}} \int _{T_e}^{T_M} {\rm d}t \langle \boldsymbol{S}_{\textit{con}}\rangle _t \nonumber \\& =\frac {1}{T_{Me}}\int _{T_e}^{T_M} {\rm d}t\frac {1}{N_p}\sum _{j=1}^{N_p}\sum _{k=1\atop {k\neq j}}^{N_p}\left (\frac {\boldsymbol{r}^j_{c,k} \boldsymbol{F}^j_{con,k}+\boldsymbol{F}^j_{con,k}\boldsymbol{r}^j_{c,k}}{2}-\frac {\boldsymbol{F}^j_{con,k}\boldsymbol{\cdot }\boldsymbol{r}^j_{c,k}\boldsymbol{I}}{3}\right ), \end{align}

where $\langle \boldsymbol{S}_{\textit{con}}\rangle _t$ is the average value of the stresslet associated with particle contact at a given time, the inner summation represents the stresslet generated by the jth particle due to mechanical contacts with neighbouring particles and its summand is the stresslet generated by the jth particle due to contact with the kth particle arising from a contact force, $\boldsymbol{F}^j_{con,k}$ , acting at the contact point $\boldsymbol{r}^j_{c,k}$ (measured relative to the centre of mass of the jth particle). At dilute concentrations, i.e. $n\ll 1$ , $\langle \boldsymbol{S}_{f}\rangle _{\textit{MI}}$ and $\langle \boldsymbol{S}_{\textit{con}} \rangle _{\textit{MI}}$ would be equivalent to $\langle \boldsymbol{S}_{f}\rangle$ and $\langle \boldsymbol{S}_{\textit{con}}\rangle$ from the PI simulations, respectively.

The transiently evolving specific viscosity is calculated as

(4.17) \begin{eqnarray} \Delta \mu (t) = n\left (\langle \boldsymbol{S}_{\textit{Hyd}}\rangle _t + \langle \boldsymbol{S}_{f}\rangle _t + \langle \boldsymbol{S}_{\textit{con}}\rangle _{t}\right ). \end{eqnarray}

Figure 13(a) illustrates the transient evolution of the specific viscosity in a suspension of SAPs or ETs of $A=40$ initially possessing an isotropic orientation distribution. Initially both SAPs and ETs exert the same stresses because most particles are in their tumbling state. The initial drop in $\Delta \mu$ for both SAPs and ETs is caused by the alignment of the particles close to the flow–vorticity plane. The ETs only temporarily align near the flow–vorticity plane, but then start tumbling and interacting with other particles through mechanical contact. In contrast, SAPs approach their equilibrium orientation $\boldsymbol{p}_S$ that is near the flow–vorticity plane closely following the transient Jeffery trajectory of isolated particles (see (1.3)), and remain near $\boldsymbol{p}_S$ , imparting much smaller shear stresses. The tilt of SAPs relative to the fluid lamellae approaches an $O(\delta )$ value in an $O(\delta ^{-1})$ time (Borker et al. Reference Borker, Stroock and Koch2024) and thereby the shear stresses also approach their steady state value of $O(\epsilon \Delta A/(A^*A^2))$ in that time as shown in figure 13(a). Understanding this transient behaviour is crucial for various applications, such as manufacturing anisotropic composite films using a coating process or fabricating a component using injection moulding. The suspension requires an $O(A\sqrt {A^*/\Delta A})$ dwell time upstream of the coating die to ensure that all particles have a sufficient time to align along the fluid lamellae. Because SAPs have weak interparticle interactions, Jeffery’s solution offers an approximate treatment of the transient change in the rheology and the orientational microstructure of the suspension, facilitating die and processing flow design.

The steady state simulations of SAPs indicated that MIs caused a minimal change in the shear stresses compared with the isolated particle result shown in figure 13(b). This primarily resulted from the tendency of SAPs to remain aligned close to the fluid lamellae even when the non-dimensional number density reaches a moderate value of $n=0.3$ . Borker et al. (Reference Borker, Stroock and Koch2024) indicated that at low particle concentrations the hydrodynamic velocity disturbance generated by SAPs was too weak to cause tumbling and at sufficiently high concentrations mechanical contact prevented the particles from tumbling causing the particles to be restricted near the flow–vorticity plane. This ensured that the value of $\langle {S}_{\textit{Hyd}}\rangle _{\textit{MI}}$ was close to $S_{\textit{Hyd}}(\boldsymbol{p}_S)$ , the shear stresslet exerted by isolated SAPs. The values of $\langle {S}_f\rangle _{\textit{MI}}$ and $\langle {S}_{\textit{con}}\rangle _{\textit{MI}}$ were found to be an order of magnitude smaller, being close to the values predicted by the PI simulations. Overall, even with the inclusion of MI effects, the shear viscosity of a suspension of SAPs was much smaller than the viscosity of ETs as shown in figure 13(b).

For a suspension of tumbling tori to have the same viscosity as an $A=40$ SAP suspension at $n=0.1$ requires the tori to have aspect ratios of O(400) based on the extrapolation of the scaling fit of $\Delta \mu =54.8n\epsilon \delta _T(1+0.9(1+3.56\epsilon )n)$ for a torus ((3.5) of Borker & Koch (Reference Borker and Koch2023)). The SAPs provide access to particle suspensions with small $\Delta \mu$ values and high degrees of anisotropy that are inaccessible using equivalent tumbling particles. The SAPs also have lower viscosity compared with other planar particles such as cylindrical discs as shown in figure 13(b). For discs, the viscosity only includes contributions from (2.2) which were evaluated using BEM simulations assuming the Leal–Hinch distribution (Leal & Hinch Reference Leal and Hinch1971). In general SAPs should have a lower viscosity compared with other planar tumbling particles. The stresslet of an isolated SAP perfectly oriented along its equilibrium orientation provided a good measure of the specific viscosity of the suspension because interparticle interactions between SAPs were weak. Therefore, $\langle S_{\textit{Hyd}}\rangle _J$ for ring-shaped SAPs with various cross-sectional shapes shown in figure 3(b) could be used to predict the specific viscosity of the suspension and one could then choose a particular SAP geometry for the desired specific viscosity. The shear stress for T-rings is approximately given by

(4.18) \begin{eqnarray} \langle S_{\textit{Hyd}}\rangle _J &=& \frac {53.5}{A^2} + \frac {125.0}{A^*}\frac {1}{A\log (8 A)} \left (\frac {A-A^*}{A}\right ), \end{eqnarray}

which captures the dominant dependence on particle aspect ratio near the critical value $A=O(A^*)$ as shown by the solid line in figure 3. Since the shear stresses for different ring-shaped SAPs are similar in magnitude, as indicated in figure 3, (4.18) also serves as a useful estimate for the particle stresslet induced by other SAP geometries in the aspect ratio range $A=O(A^*)$ .

Figure 13.

Shear stress from MI simulations. (a) Transient evolution of the specific viscosity in a suspension of SAPs and ETs, with particles initially oriented isotropically. The specific viscosity of SAP suspensions decreases rapidly as the particles align with their equilibrium orientation. (b) Specific viscosity as a function of the non-dimensional number density $n$ for SAPs, ETs and discs. For SAPs, the specific viscosity is obtained from isolated particle calculations using (2.2), and from PI and MI simulations. For ETs, results from PI and MI simulation are presented. For discs, only the isolated particle values are included and calculated using (2.2). The errorbars and shaded regions represent 95 % confidence intervals.

4.4. First and second normal stress difference

As discussed in § 2, SAP suspensions exhibit non-zero normal stress differences even at dilute concentrations due to the asymmetry in their orientation relative to the gradient–vorticity plane. In this subsection, we include the effects of interparticle interactions on normal stress differences and demonstrate its influence on the suspension flow. The first and second normal stress difference in the suspension calculated from PI simulations are given by

(4.19) \begin{eqnarray} N_{i} &=& \langle N_i\rangle _J + \langle N_i\rangle = \langle N_i\rangle _J + n\bigl ( \delta _{i1} (\hat {\boldsymbol{e}}_1\hat {\boldsymbol{e}}_1- \hat {\boldsymbol{e}}_2\hat {\boldsymbol{e}}_2)\nonumber \\ &+&\delta _{i2}(\hat {\boldsymbol{e}}_2\hat {\boldsymbol{e}}_2 - \hat {\boldsymbol{e}}_3\hat {\boldsymbol{e}}_3)\bigr ):\left ( \langle \boldsymbol{S}_{p}\rangle + \langle \boldsymbol{S}_f\rangle + \langle \boldsymbol{S}_{\textit{con}} \rangle \right ), \end{eqnarray}

and from MI simulations is

(4.20) \begin{align} N_{i} =&\langle N_i\rangle _{\textit{MI}} = n\bigl ( \delta _{i1} (\hat {\boldsymbol{e}}_1\hat {\boldsymbol{e}}_1- \hat {\boldsymbol{e}}_2\hat {\boldsymbol{e}}_2)\nonumber \\ &\quad+\delta _{i2}(\hat {\boldsymbol{e}}_2\hat {\boldsymbol{e}}_2 - \hat {\boldsymbol{e}}_3\hat {\boldsymbol{e}}_3)\bigr ):\left ( \langle \boldsymbol{S}_{\textit{Hyd}}\rangle _{\textit{MI}}+ \langle \boldsymbol{S}_f\rangle _{\textit{MI}}+ \langle \boldsymbol{S}_{\textit{con}} \rangle _{\textit{MI}} \right ), \end{align}

where $i\in \{1,2\}$ . The $O(n^2)$ contributions to the normal stress differences from PIs, i.e. $\langle N_i\rangle$ , were notably smaller for SAPs compared with ETs because of their weaker interparticle interactions. Figure 14(a,b) illustrates the variation of $\langle N_1\rangle /n^2$ and $\langle N_2\rangle /n^2$ with the particle aspect ratio using PI simulations. Although SAPs exhibit non-zero normal stresses at $O(n)$ , our simulations indicated that their interactions induced minimal alterations relative to the isolated particle values described in figure 5(a). In contrast, normal stresses in a suspension of ETs were generated from colliding trajectories (Borker & Koch Reference Borker and Koch2023). Colliding trajectories for tumbling rings induced asymmetry in the particle’s orientational trajectory about the gradient–vorticity and flow–vorticity planes and generated a much larger $O(\epsilon )$ stresses during the interaction.

Figure 14(c,d) describe the variation of the first and second normal stress difference in a suspension of SAPs and ETs as a function of the non-dimensional number density, $n$ . These results confirm that the magnitudes of normal stresses generated by SAPs are numerically higher than those generated by ETs at the same particle concentrations, in the range $0\lt n\lt 0.3$ . Moreover, the periodic box simulations (symbols) for SAPs corroborated the finding that interparticle interactions contributed little to the normal stresses and the major contribution arose from asymmetry of the equilibrium orientation about the gradient–vorticity plane. This indicates that the normal stress difference in SAP suspensions can be accurately estimated using isolated-particle BEM calculations (Borker et al. Reference Borker, Stroock and Koch2018), even at moderate concentrations. For T-ring shaped SAPs, the normal stress differences are well described by the semiempirical relations

(4.21) \begin{align} \langle N_{1}\rangle _J = n\frac {80.2}{A\log (8A)} \sqrt {\frac {A-A^*}{A^*}}, \\[-28pt] \nonumber \end{align}

(4.22) \begin{align} \langle N_{2}\rangle _J = -0.22 \langle N_{1}\rangle _J. \\[0pt] \nonumber \end{align}

This functional form also applies to other ring-shaped SAPs, with the prefactor in (4.21) lying in the range $60{-}120$ .

Figure 14.

Contributions to the normal stresses from interparticle interactions. (a) The PI contributions to the first ( $\langle N_1\rangle$ ) and (b) second ( $\langle N_2\rangle$ ) normal stress difference as a function of the particle aspect ratio for T-rings ( $\square$ ) and tori (o). The solid line represents extrapolation of the scaling fit for tori from Borker & Koch (Reference Borker and Koch2023). (c) The first ( $N_{1}$ ) and (d) second ( $N_{2}$ ) normal stress difference for a suspension of $A=40$ SAP (–) and ET (:) using PI simulations (lines) and using MI simulations (symbols). Shaded regions and error bars are 95 % confidence intervals from the corresponding simulation values.

The magnitude of the normal stress difference provides a direct measure of the non-Newtonian nature of a particle suspension. When normal stress differences are negligible, the rheology of particle suspensions can be adequately described by an effective viscosity. The extent of the non-Newtonian behaviour is commonly quantified by the ratio of first normal stress difference to the shear stress in the suspension, defined as $\alpha _i = N_i/(1 +\sigma _{p,12})$ (Singh & Nott Reference Singh and Nott2003; Maklad & Poole Reference Maklad and Poole2021). In SAP suspensions, $\alpha _1$ is positive, while $\alpha _2$ is negative and approximately one-fifth the magnitude of $\alpha _1$ . This trend is qualitatively similar to that observed in suspensions of other high-aspect-ratio, tumbling particles such as fibres (Petrich, Koch & Cohen Reference Petrich, Koch and Cohen2000; Snook et al. Reference Snook, Davidson, Butler, Pouliquen and Guazzelli2014; Bounoua et al. Reference Bounoua, Kuzhir and Lemaire2016) and tori (Borker & Koch Reference Borker and Koch2023). At sufficiently low particle volume fractions, $\varPhi$ , SAP suspensions generate larger normal stress differences than suspensions of other rigid particles such as spheres (Zarraga et al. Reference Zarraga, Hill, Leighton and David2000; Singh & Nott Reference Singh and Nott2003), fibres (Snook et al. Reference Snook, Davidson, Butler, Pouliquen and Guazzelli2014; Bounoua et al. Reference Bounoua, Kuzhir and Lemaire2016) or discs (Bertevas, Fan & Tanner Reference Bertevas, Fan and Tanner2010). This behaviour is illustrated in figure 15(a), which describes the variation of $\alpha _1$ with the particle volume fraction. For visual clarity, the data for sphere suspensions are plotted as $-\alpha _1$ to facilitate comparison. Experimental measurements indicate that, at sufficiently low volume fractions, suspensions of spheres or fibres generate negligible normal stresses due to dominant HIs that preserve microstructural symmetry. The dashed and dotted lines extrapolate the PI simulation results to higher concentrations ( $0.3\lt n\lt 2$ ) for $A=40$ SAP and ET suspensions, respectively. These trends suggest that SAP suspensions can sustain higher $\alpha _1$ values than other high-aspect-ratio particles. Only the experimental results for a suspension of $A=50$ fibres by Petrich et al. (Reference Petrich, Koch and Cohen2000) approach the values predicted for SAPs at similar particle volume fractions. However, their analysis assumed $N_2=0$ , implying that their values depicted in figure 15(a) represent $(N_1-N_2)/(1 + \sigma _{p,12})$ . Since more recent measurements have shown that $N_2$ is typically negative and comparable in magnitude to $N_1$ (Snook et al. Reference Snook, Davidson, Butler, Pouliquen and Guazzelli2014; Bounoua et al. Reference Bounoua, Kuzhir and Lemaire2016), their value of $\alpha _1$ is likely overestimated compared with SAP suspensions. Figure 15(b) also shows the ratio of the second and the first normal stress difference, highlighting that the magnitudes of $N_1$ and $N_2$ are comparable across most particle suspensions. Together, these results demonstrate that SAP suspensions exhibit non-Newtonian effects of similar magnitude as other particle suspensions, but at much lower particle volume fractions.

Figure 15.

(a) Ratio of the first normal stress difference to the shear stress, and (b) ratio of the second to the first normal stress difference, as functions of particle volume fraction for suspensions of SAPs, tori, fibres, discs, and spheres. For SAPs, data are shown from MIs () and PIs (). The dashed (–) and dotted (:) lines are extrapolations of the PI results for $A=40$ SAPs and tori (Borker & Koch Reference Borker and Koch2023), respectively, with the shaded regions representing 95 % confidence intervals. Other datasets include: $A=40$ tori (Borker & Koch Reference Borker and Koch2023), $A=18, 33$ fibres (Bounoua et al. Reference Bounoua, Kuzhir and Lemaire2016), $A=17, 32$ fibres (Snook et al. Reference Snook, Davidson, Butler, Pouliquen and Guazzelli2014), $A=50$ fibres (Petrich et al. Reference Petrich, Koch and Cohen2000), $A=5$ discs (Bertevas et al. Reference Bertevas, Fan and Tanner2010) and spheres ( $\diamond$ = Singh & Nott (Reference Singh and Nott2003)). The data for spheres represents $-\alpha _1$ in (a) and $N_2/N_1$ in (b) for ease of visual representation. Additional experimental data extending beyond the plotted range can be found in the respective references.

Non-zero normal stress differences in particle suspensions can lead to various phenomena, including particle migration (Morris & Boulay Reference Morris and Boulay1999; Ramachandran & Leighton Jr Reference Ramachandran, Leighton and David2008), the Weissenberg effect (rod climbing or dipping phenomenon) (Weissenberg Reference Weissenberg1947; Nunez et al. Reference Nunez, Ribeiro, Arney, Feng and Joseph1994; Boyer, Pouliquen & Guazzelli Reference Boyer, Pouliquen and Guazzelli2011; More et al. Reference More, Patterson, Pashkovski and McKinley2023), secondary flows in curved or non-circular pipes (Barnes Reference Barnes2003; Maklad & Poole Reference Maklad and Poole2021), interfacial instabilities involving flow of multiple particle suspensions (Brady & Carpen Reference Brady and Carpen2002) and die swelling (Tanner Reference Tanner1970; Mezi et al. Reference Mezi, Ausias, Grohens and Férec2019). While these effects are prominent in polymeric liquids, they can also emerge in particle suspensions at sufficiently large particle concentrations (Morris & Boulay Reference Morris and Boulay1999; Boyer et al. Reference Boyer, Pouliquen and Guazzelli2011). Given that SAP suspensions generate significant normal stresses, many of these effects should be observable even at moderate particle concentrations of $n=O(1)$ . As an example, we examine the Weissenberg effect by estimating the meniscus rise when a rod rotates in a large pool of SAP suspension.

The Weissenberg effect is most notable in polymeric liquids (Barnes, Hutton & Walters Reference Barnes, Hutton and Walters1989; More et al. Reference More, Patterson, Pashkovski and McKinley2023), where a positive first normal stress difference acts as a hoop stress around the rotating rod, driving fluid radially inward and causing the meniscus to rise until it reaches a steady height profile (Barnes et al. Reference Barnes, Hutton and Walters1989). Dense suspensions of spheres, which exert negative first and second normal stress differences, generate a rod-dipping effect instead (Zarraga et al. Reference Zarraga, Hill, Leighton and David2000; Boyer et al. Reference Boyer, Pouliquen and Guazzelli2011), but only at relatively high concentrations ( $\varPhi \gt 0.22$ ). Fibre suspensions, with a positive first and negative second normal stress difference, have shown minimal meniscus deformations under similar conditions (Snook et al. Reference Snook, Davidson, Butler, Pouliquen and Guazzelli2014). Our simulations suggest that SAP suspensions produce sufficiently large normal stress differences for rod climbing to be measurable at moderately large wall shear rates of $100$ $\textrm {s}^{-1}$ , which is within the range of recent experiments (Boyer et al. Reference Boyer, Pouliquen and Guazzelli2011; Snook et al. Reference Snook, Davidson, Butler, Pouliquen and Guazzelli2014).

Figure 16.

Weissenberg effect in SAP suspensions. (a) The meniscus height of SAP suspensions relative to the height far away from the rod, $h$ , as a function of the radial distance from the rod surface for T-ring shaped SAPs of different aspect ratios (solid lines) and $A=40$ tori. Here, $\mu_{\kern-1pt f} =2$ Pa s, $\varOmega = 100$ $\textrm {s}^{-1}$ , $\rho =1260$ kg $\textrm {m}^{-3}$ , $n=0.2$ and $\textit{Pe}\gg \delta ^{-3}$ . Here SAPs with $A=28, 102$ only include the isolated particle contributions described in figure 5(a). Here SAP and ET with $A=40$ additionally included the PI contributions and the shaded region represents 95 % confidence interval. (b) Variation of the height of the suspension surface at the rod, $h_0$ , in millimetres as a function of the aspect ratio, $A$ , and non-dimensional number density $n$ .

Following the derivation in Zarraga et al. (Reference Zarraga, Hill, Leighton and David2000) in the absence of surface tension and inertial forces, the height profile of the suspension meniscus, $h$ , relative to the height far from the rod is given by

(4.23) \begin{eqnarray} h(r) = \left (\frac {{N}_{1}}{2} + {N}_{2}\right ) \frac {2\mu_{\kern-1pt f}\varOmega }{\rho g} \frac {b^2}{r^2}, \end{eqnarray}

where $r$ is the radial distance from the axis of the rotating rod, ${N}_{1}$ and ${N}_{2}$ are the non-dimensional first and second normal stress differences described earlier, $\varOmega$ is the rotational speed of the rod, $\rho$ is the suspension density, $b$ is the radius of the rod and $g=9.81$ m $\textrm {s}^{-2}$ is the gravitational acceleration. Figure 16(a) presents the profile of the meniscus for SAP suspensions at $n=0.2$ for various aspect ratios, using glycerine as the background fluid ( $\mu_{\kern-1pt f} = 2$ Pa s, $\rho = 1260$ kg $\textrm {m}^{-3}$ ), and a shear rate $\gamma = 100$ $\textrm {s}^{-1}$ . The profile for an $A=40$ ET is also depicted (dashed line) to demonstrate that the rise in the meniscus height for a suspension of tumblers is significantly lower at the same $n$ . The meniscus height at the rod surface, termed $h_{0}$ , is then

(4.24) \begin{eqnarray} h_{0} = \left (\frac {{N}_{1}}{2} + {N}_{2}\right ) \frac {2\mu_{\kern-1pt f}\varOmega }{\rho g}. \end{eqnarray}

Because SAPs experience negligible interparticle interactions, as described earlier, the isolated particle results from § 2 can be used, giving ${N}_{1} = n\xi _3 \sin (4\delta )/4$ , and ${N}_{2}=-0.22{N}_{1}$ . Substituting this in (4.23), yields the following meniscus profile for SAP suspensions:

(4.25) \begin{eqnarray} h(r) = \frac {0.56 N_1\mu_{\kern-1pt f}\varOmega }{\rho g} \frac {b^2}{r^2}= \frac {0.14n\mu_{\kern-1pt f}\varOmega \xi _3\sin (4\delta )}{\rho g} \frac {b^2}{r^2}. \end{eqnarray}

The values of $N_1$ for various cross-sectional shapes are depicted in figure 5(b). The values of $N_1$ , (or $\xi _3$ and $\delta$ ) for other cross-sectional shapes can be estimated from SBT (Borker & Koch Reference Borker and Koch2019) or BEM (Borker et al. Reference Borker, Stroock and Koch2018, Reference Borker, Stroock and Koch2024) computations. Figure 16(b) shows the variation of $h_0$ with the aspect ratio, $A$ and non-dimensional number density $n$ for T-rings. At low number densities, $h_0$ is typically larger for SAPs than tumbling T-rings, with a peak value occurring near $A\approx 55(\approx 2 A^*)$ . For tumbling T-rings with $A\ll A^*$ , $h_0$ also rises with increasing $n$ due to increased frequency of interparticle interactions. We observe that at the same aspect ratio T-rings also produced larger rod climbing heights than ETs for $n=O(1)$ . The predicted rod climb heights for a suspension of SAPs in glycerine at angular velocities of $O(100)$ $\textrm {s}^{-1}$ and $n=O(1)$ are comparable to the lower range of values measured for polymer solutions (More et al. Reference More, Patterson, Pashkovski and McKinley2023), suggesting that rod climbing should be observable under practical laboratory conditions (Boyer et al. Reference Boyer, Pouliquen and Guazzelli2011; More et al. Reference More, Patterson, Pashkovski and McKinley2023). This makes experimental validation of the pronounced non-Newtonian behaviour of SAP suspensions feasible.