Stochastic resonance of rotating particles in turbulence
Ziqi Wang, Xander M. de Wit, Roberto Benzi, Chunlai Wu, Rudie P. J. Kunnen, Herman J. H. Clercx, Federico Toschi

TL;DR
Researchers found that magnetic particles can detect and influence small-scale turbulence through a phenomenon called stochastic resonance, opening new ways to study and control turbulent flows.
Contribution
The study introduces the first experimental and numerical observation of stochastic resonance in turbulent particle dynamics and a novel symmetry-breaking mechanism.
Findings
Magnetic particles exhibit stochastic resonance in turbulence, enhancing rotational response to external magnetic fields.
A symmetry-breaking mechanism allows net particle rotation in zero-mean vorticity turbulence using oscillating magnetic fields.
The findings suggest a magnetic resonance-based method for measuring turbulent vorticity in optically inaccessible conditions.
Abstract
The chaotic dynamics of small-scale vorticity plays a key role in understanding and controlling turbulence, with direct implications for energy transfer, mixing, and coherent structure evolution. However, measuring or controlling its dynamics remains a major conceptual and experimental challenge due to its transient and chaotic nature. Here we use a combination of experiments, theory, and simulations to show that small magnetic particles of different densities, exploring flow regions of distinct vorticity statistics, can act as effective probes for measuring and forcing turbulence at its smallest scale. The interplay between the magnetic torque, from an externally controllable magnetic field, and hydrodynamic stresses, from small-scale turbulent vorticity, uncovers an extremely rich phenomenology. Notably, we present the first experimental and numerical observation of stochastic…
Genes, proteins, chemicals, diseases, species, mutations and cell lines named across the full text — each resolved to its canonical identifier and authoritative record.
Click any figure to enlarge with its caption.
Figure 1
Figure 2
Figure 3
Figure 4
Figure 5- —501100003246Nederlandse Organisatie voor Wetenschappelijk Onderzoek (Netherlands Organisation for Scientific Research)
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Taxonomy
Topicsstochastic dynamics and bifurcation · Particle Dynamics in Fluid Flows · Micro and Nano Robotics
Introduction
Turbulent flows, characterized by chaotic multi-scale fluctuations, are central to a wide range of natural and industrial processes^1,2^. A key challenge in turbulence research is the direct measurement and control of the vorticity, which is primarily concentrated in small-scale vortex filaments^3,4^. Due to their transient nature and rapid chaotic evolution, capturing the rotational dynamics of these filaments remains a major conceptual and experimental challenge. Several methods have been developed to infer small-scale turbulence characteristics, including particle tracking methods (e.g., tracer particles, deformable particles, and pattern-coated particles), as well as thermal and optical anemometry^5–7^. However, these approaches often suffer from inherent limitations such as insufficient spatial and temporal resolution, signal occlusion in optically dense environments, and difficulties in directly resolving vorticity at the smallest turbulence scales. As a result, achieving real-time, high-fidelity measurements of rotational motion in turbulent flows remains an ongoing challenge, particularly in complex and optically inaccessible flows.
A promising alternative is the use of magnetic particles, which respond to both stochastic vorticity fluctuations of turbulence (hydrodynamic torques^8^) and deterministic torques imposed by external magnetic fields. The dynamics of buoyant magnetic particles, in particular, are closely linked to the local vorticity. Particles with different densities exhibit distinct rotational behaviors as they explore different regions of the flow, leading to preferential concentration: light (heavy) particles tend to accumulate in high- (low-) vorticity regions^7,9,10^. These particles offer a unique dual role: they can serve as candidates for probing rotational fluctuations, while also enabling active flow manipulation through external forcing. This ability to modulate particle rotation and measure its response in real time provides a new pathway for investigating turbulence at small scales and suggests an elegant strategy for measuring turbulent vorticity. However, in order to unlock these novel directions, an accurate understanding of the complex interplay between hydrodynamic and magnetic torques on particle dynamics is required.
Previous studies have extensively explored magnetic particles in quiescent fluids under rotating magnetic fields, revealing rich phenomena such as self-assembly into chain-like structures^11,12^, synchronization-selected structures^13^, and more complex arrays^14^, enhanced mixing^11,12,15^, turbulence design^16–18^ and phase locking, where particles follow the rotating magnetic field with a constant phase lag with respect to the magnetic field when the rotating magnetic angular velocity ωH is below a threshold ωcr^19^. Investigations have also extended to self-propelling active magnetic particles (e.g., magnetotactic bacteria)^20,21^, magnetic Janus particles^13,22^, particles with different aspect ratios^23,24^, colloidal spinners and active swimmers^25–27^.
The competition and interplay between effects of turbulence and magnetic fields may lead to complex particle dynamics, yet to be understood, presenting a fundamental challenge in predicting and controlling particle behavior in turbulent environments.
Here, we explore this question by investigating the rotational dynamics of magnetic particles in turbulence under a rotating magnetic field. While deterministic synchronized and oscillation (back-and-forth) rotational regimes have been observed in zero-noise systems (e.g., compasses or rolling particles in quiescent fluids^28–30^), our work extends this understanding into the complex environment of fully developed turbulence. We systematically explore the complex interplay between deterministic external magnetic forcing and turbulent hydrodynamic torques, revealing a significantly richer phase space of particle rotational dynamics. We show that the effect of turbulence (vorticity) acts as an effective noise (Supplementary Information) and reveal the first signature of the emergence of stochastic resonance (SR)—a counterintuitive phenomenon where internal fluctuations non-linearly cooperate with external forcing to amplify the system response^31–39^—of particles in turbulence, which can be used as the mechanism of magnetic particles acting as effective vorticity probes for a novel type of turbulence “microscope”. While the concept of SR, initially observed in bistable systems subjected to noise and weak periodic signals, has since been extended to a wide range of physical and biological systems, as well as signal processing and computer science^40–45^, its role in turbulence-driven rotational dynamics remains mostly unexplored.
To identify potential SR signatures, we set out to leverage turbulent vorticity fluctuations as a source of noise and examine the response of the magnetic particle dynamics to the rotating magnetic field, by means of experiments, simulations, and theoretical modeling. Our experimental platform (Fig. 1a and e) involves a dilute collection of approximately spherical light magnetic particles (produced by Styrofoam core coated with magnetic paint, Fig. 1b, with a mean density of 20% with respect to water and with a mean diameter of 0.7 mm, see “Methods” for details) which are immersed in a Von Kármán-type turbulent flow (producing homogeneous and isotropic turbulence (HIT) in the center region of the cell^46–49^) in an octagonal water-filled chamber. Magnetic particles are designed to be light and of a size comparable to the Kolmogorov scale, allowing them to preferentially enter and interact with vortex filaments. The flow is generated by two counter-rotating bladed disks and the chamber is placed in a system of Helmholtz coils, which is programmed to generate a uniform magnetic field in a restricted measurement volume, H = H h, rotating in the xoy-plane with angular velocity vector \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\boldsymbol{\omega }}}}}_{H}={\omega }_{H}\hat{{{{\bf{z}}}}}$$\end{document} where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hat{{{{\bf{z}}}}}$$\end{document} is the unit vector along the z-axis (“Methods” and Supplementary Information). As the induced magnetic dipole, m, of a particle tends to align with the magnetic field direction, h, the particle experiences a magnetic torque balanced by the rotational drag and spins around an axis in the particle body frame parallel to \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hat{{{{\bf{z}}}}}$$\end{document} along its trajectory (as shown by typical stroboscopic time-lapse trajectories in Fig. 1c, d). The rotational motion is measured by tracking particle surface patterns (“Methods”).Fig. 1. Rotational dynamics of magnetic particles in a rotating magnetic field: “phase-locked” versus “back-and-forth” regimes.a Magnetic particles in turbulence under a rotating magnetic field: The experimental setup consists of a turbulence generator, a magnetic field generator, and magnetic particles. A Von Kármán-type turbulent flow is generated in an octagonal water-filled container with an internal diameter of 2R = 150 mm and a height of 220 mm, driven by two counter-rotating bladed disks. A uniform rotating magnetic field H = Hh, with angular frequency ωH, is generated using a system of two pairs of perpendicular Helmholtz coils. b Magnetic particles consist of a Styrofoam core coated with magnetic paint. Surface patterns enable tracking rotational motion, as illustrated in the stroboscopic time-lapse trajectories (c, d). e Photograph of the setup. f–i In a quiescent fluid, particles exhibit two rotational regimes: “phase-locked” (low ωH) and “back-and-forth” (high ωH). f (Supplementary Movie 1) and h (Supplementary Movie 2) show rotational trajectories within one rotation period, viewed in a plane perpendicular to the rotation plane of the magnetic field. The particle initial preferred magnetization direction n0 (dashed arrows) and instantaneous orientation n (solid arrows) are marked over time. During rotation within a quiescent fluid, n and h remain within the same plane^50^. The corresponding normalized particle angular velocity, ωp,i/ωH (i = x, y, z), is shown in (g and i). Detailed information about the experiments and simulations can be found in “Methods”.
Results
Theoretical representation
The particle rotation is characterized by a preferred magnetization direction, n (Fig. 2a), whose trajectory can be tracked in turbulence under the influence of a rotating magnetic field. By balancing the magnetic torque and drag torque induced by turbulence^19,50^ (“Methods” and Supplementary Information), the equations governing the particle rotational motion can be obtained, in the overdamped limit, as
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\frac{d{{{\bf{n}}}}}{dt}={\omega }_{{{{\rm{a}}}}}({{{\bf{n}}}}\cdot {{{\bf{h}}}})\left({{{\bf{h}}}}-({{{\bf{n}}}}\cdot {{{\bf{h}}}}){{{\bf{n}}}}\right)+{{{{\boldsymbol{\omega }}}}}_{{{{\rm{f}}}}}\times {{{\bf{n}}}},$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\frac{d\left({\xi }_{{{{\rm{r}}}}}{{{{\boldsymbol{\omega }}}}}_{{{\rm{p}}},\parallel }\right)}{dt}=-{\xi }_{{{{\rm{r}}}}}\left({{{{\boldsymbol{\omega }}}}}_{{{\rm{p}}},\parallel }-\left({{{{\boldsymbol{\omega }}}}}_{{{{\rm{f}}}}}\cdot {{{\bf{n}}}}\right){{{\bf{n}}}}\right),$$\end{document}where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\omega }_{{{{\rm{a}}}}}=\mu {{{{\mathcal{V}}}}}_{{{{\rm{p}}}}}{H}^{2}({\chi }_{\parallel }-{\chi }_{\perp })/{\xi }_{{{{\rm{r}}}}}$$\end{document} is the magnetic field strength normalized by the rotational drag coefficient ξr, with μ, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\mathcal{V}}}}}_{{{{\rm{p}}}}}$$\end{document} , and χ∥,⊥ the particle permeability, volume, and anisotropic magnetic susceptibilities, respectively^19,50,51^ (“Methods”). The term ωp,∥ represents the particle spinning rate. The perpendicular rotational component, ωp,⊥, can be obtained by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\boldsymbol{\omega }}}}}_{{{\rm{p}}},\perp }={{{\bf{n}}}}\times \dot{{{{\bf{n}}}}}$$\end{document} which represents the particle tumbling rate. The turbulent vorticity experienced by the particle is ωf, which is defined as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\boldsymbol{\omega }}}}}_{{{{\rm{f}}}}}=\frac{1}{2}{[{{{\boldsymbol{\nabla }}}}\times {{{{\bf{u}}}}}_{{{{\rm{f}}}}}]}_{{{{\rm{O}}}}}$$\end{document} with the subscript O implying evaluation at the center of mass of the particle. The total particle angular velocity is thus ωp = ωp,∥ + ωp,⊥.Fig. 2. Experimental validation of the theoretical model for particle rotational dynamics with turbulence.a Schematic of a magnetic particle in turbulence under a rotating magnetic field. The particle with preferred magnetization direction n is subjected to a magnetic field h in the xoy-plane, rotating around the z-axis with frequency ωH, and experiences turbulent vorticity ωf in homogeneous isotropic turbulence. The phase lag β (defined as the angle between the projection of n, i.e., nxoy, and h) and the rotation angle ψ (angle of nxoy relative to the x-axis) describe the particle dynamics. The dynamical evolution of β is governed by the interplay of the magnetic field rotation drag, magnetic torque, and turbulence torque. b–d Comparison of the probability density distribution function (PDF) of particle angular velocity in experimental (circles) and numerical (squares) studies. The experimental results (circles) are shown for varying magnetic field strengths with fixed turbulence intensity. b weak (1.2 mT, Supplementary Movie 3), c intermediate (1.4 mT, Supplementary Movie 4), and d strong (1.6 mT, Supplementary Movie 5), with constant rotational frequency of the magnetic field (parameter settings can be found in the “Methods” and Supplementary Information).
In a noiseless system, i.e., ωf = 0, the critical frequency is ωcr = ωa/2 (details in Supplementary Information). For ωH ≤ ωcr, the particle is synchronized with the rotating magnetic field ("phase-locked”), i.e., ωp,z = ωH. Otherwise, when ωH > ωcr, the particle undergoes back-and-forth motion ("back-and-forth”) when the maximum magnetic torque cannot balance the response drag torque: the phase lag angle β between nxoy (the projection of n on xoy-plane) and h changes periodically over time^20^. These two regimes can be observed experimentally (Fig. 1f–i). In the “phase-locked” regime (Supplementary Movie 1), the particle rotates synchronously with the magnetic field, as evidenced by the normalized angular velocity, ωp,z/ωH ≈ 1, in Fig. 1g. In the “back-and-forth” regime (Supplementary Movie 2), the particle oscillates periodically (Fig. 1i).
As illustrated in the schematic diagram of Fig. 2a, the dynamics of the particle preferred magnetization direction, n, is fundamentally characterized by the phase lag angle, β, with respect to h on the *xoy-*plane. The governing equations, Eqs. (1), can be simplified and explicitly written as (derivations in “Methods”)
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{\beta }={\omega }_{H}-\frac{{\omega }_{{{{\rm{a}}}}}}{2}\sin 2\beta -{\omega }_{{{\rm{f}}},z}.$$\end{document}Geometrically, we have β = ωHt–ψ. The phase lag angle β represents the particle rotational displacement in the co-rotating frame of the magnetic field, while ψ is the rotation angle of the projection n on the xoy-plane with respect to the x-axis, which describes the particle rotational displacement in the laboratory frame. The turbulent vorticity signal, ωf, represents the vorticity experienced by the particle center along its trajectory. This simplified representation, Eq. (2), separates different contributions and clearly shows that the rotational motion of the particle results from the interplay between the magnetic field rotational drag, the magnetic torque, and the turbulence drag torque (Fig. 2a).
Switching on turbulence
Now we wonder how this picture will change when turbulence is taken into consideration. The turbulence induces fluctuations on the small particles of a typical magnitude ωη = 1/τη, the Kolmogorov characteristic frequency. Experimentally, we turn on the turbulence and vary the intensity of the magnetic field. By statistically analyzing the probability density function (PDF) of the particle angular velocity ωp, a dominant peak is observed in the z-direction at approximately ωp,z ≈ ωH (Fig. 2b–d, purple circles), while the angular velocities in the x and y directions are similar to each other but significantly different from that in the z-direction (Fig. 2b–d, yellow and green circles).
To understand the influence of turbulence on the particle dynamics further, we carry out a systematic numerical investigation. To model the motion of light particles in HIT, we first perform direct numerical simulations using a pseudo-spectral method^10,52,53^ to obtain ωf (see “Methods” for details). The resulting vorticity data is then incorporated into the theoretical model, Eqs. (1), which is integrated using a third-order Runge–Kutta method.
The simulated particle rotation is validated against the experimental results (Fig. 2b–d and Supplementary Movies 3–5, respectively). Notably, the theoretical model (squares) accurately captures the essential statistical physics (peak positions) of the particle stochastic rotational dynamics in the experiments (circles). The slight discrepancies in the PDF shapes arise from polydispersity in particle shape and magnetic properties in experiments which influences how individual particles respond to an external magnetic field, whereas the simulations assume idealized spherical particles with the same magnetic characteristics for each particle (“Methods” and Supplementary Information). Further comparisons between simulation and experimental results all demonstrate good qualitative agreement (Supplementary Information, Supplementary Movies 9 and 10).
Additionally, to gain deeper insight into the stochastic nature of particle rotation in turbulence, we analyze the waiting time, τ—the duration between successive transitions of β between its metastable states—to characterize the nature of turbulence-induced noise in the absence of external forcing (i.e., ωH = 0). The Poisson-like distribution of τ confirms that, within the investigated parameter space, stochastic transitions induced by turbulent (vorticity) fluctuations can be effectively modeled using white noise with an equivalent intensity (Supplementary Information).
Three distinct regimes of particle dynamics
With the theoretical model capturing the core dynamics of this complex system, and with simulations allowing exploration of a broader parameter space and finer details of particle rotation, we systematically investigate how turbulence-induced vorticity fluctuations and the deterministic magnetic forcing shape particle rotational dynamics.
Beyond the previously identified “phase-locked” and “back-and-forth” regimes, we uncover a third regime: the turbulence-dominated regime, observed at small ωa w.r.t. ωη, when the effect from the magnetic field is negligible. Here, the particle rotation is primarily driven by turbulence, and its normalized angular velocity, ωp/ωH, fluctuates randomly around zero with no net angular drift (Fig. 3b), as evidenced by the random orientation distribution of the particle orientation, n (Fig. 3d and Supplementary Movie 6). Consequently, the phase lag β and its derivative are also randomized (Fig. 3c). In the “phase-locked” regime (ωa ≫ ωη, ωH < ωcr), the particle phase lag, β, remains constant, meaning \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{\beta }=0$$\end{document} (Fig. 3g). The particle angular velocity along the magnetic field rotation axis (z-axis) is locked to ωH (Fig. 3f, green line). The particle orientation is strongly confined, exhibiting 2D behavior (Fig. 3h and Supplementary Movie 7). In the “back-and-forth” regime (ωa ≫ ωη, ωH > ωcr), the particle undergoes periodic acceleration and deceleration (Fig. 3j, k, inset), forming distinct looping patterns in its orientation trajectory (Fig. 3l and Supplementary Movie 8). These regimes are summarized in a phase diagram, characterized by the normalized variance of the particle angular velocity in the z-direction, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle {\omega }_{\,{{\rm{p}}}\,,z}^{2}-{\langle {\omega }_{{{\rm{p}}},z}\rangle }^{2}\rangle /{\omega }_{\eta }^{2}$$\end{document} , where 〈 ⋅ 〉 represents a time average. The noiseless critical frequency ωcr = ωa/2 is marked by a black dot-dashed line. Unlike in a quiescent fluid, the phase-locked regime (Fig. 3m, purple region) is bounded by the blue line, predicted by the scaling argument based on the simplified governing equation of Eq. (2): ωa = 2ωH + ωη (detailed reasoning can be found in “Methods”). This shift in the subcritical boundary, relative to the noiseless case (ωa = 2ωH), accounts for the influence of turbulence fluctuations, which have a standard deviation on the order of ωη.Fig. 3. Three distinct regimes of magnetic particle rotation dynamics in turbulence under a rotating magnetic field.a Turbulence-dominated regime (ωa ≪ ωη). b The normalized particle angular velocity, ωp/ωH, exhibits turbulent fluctuations with a zero-average and the phase lag derivative, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{\beta }$$\end{document} is randomized (c). d The evolution of the tip of preferred magnetization direction of the particle, n, is visualized in space, with color indicating time progression (Supplementary Movie 6). The random distribution of orientations confirms the turbulent-dominated nature. e “Phase-locked” regime (ωa ≫ ωη, ωH < ωcr): The angular velocity of the particle along the magnetic field rotation axis (z-axis) is locked to ωH (f, green line), with zero net drift for other components (f, blue and orange lines). The phase lag derivative \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{\beta }$$\end{document} is effectively zero (g) and the particle orientation follows a 2D trajectory (h, Supplementary Movie 7). i “Back-and-forth” regime (ωa ≫ ωη, ωH > ωcr): The angular velocity of the particle along the z-axis exhibits periodic acceleration and deceleration oscillations (j, inset), and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{\beta }$$\end{document} varies similarly (k, inset). The orientation vector oscillates periodically (back-and-forth looping in l, Supplementary Movie 8). m, Phase diagram of particle rotation regime, colored by the normalized variance σ(ωp,z) of the particle angular velocity along the z-direction, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sigma ({\omega }_{{{\rm{p}}},z})=\langle {\omega }_{\,{{\rm{p}}}\,,z}^{2}-{\langle {\omega }_{{{\rm{p}}},z}\rangle }^{2}\rangle /{\omega }_{\eta }^{2}$$\end{document} . Black dot-dashed line: the critical frequency, ωcr = ωa/2. Green dashed line: the condition where the turbulence fluctuation intensity is balanced by the magnetic strength, i.e., ωa = ωη. The phase-locked regime is bounded with the blue thick line, predicted by the scaling analysis, i.e., ωa = 2ωH + ωη. Symbols (square, star, and circle): representative cases from the three distinct regimes shown in earlier panels. Triangle markers: experimental settings from Fig. 2b–d. The diamond symbols denote the experimental results at ωH/ωa = 0.12 ± 0.01, providing evidence of stochastic resonance, which will be discussed further in Fig. 4a and c.
Stochastic resonance
To determine whether SR occurs, we focus on the β and its behavior across different parameter settings. We summarize the time-averaged value of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }\rangle /{\omega }_{\eta }$$\end{document} in one phase diagram (Fig. 4a), with ωη/ωa (vertical axis) quantifying the relative turbulent fluctuation intensity (noise) to the applied magnetic field strength, and ωH/ωa (horizontal axis) distinguishing noiseless subcritical, critical (i.e., ωH/ωa = 1/2, black dot-dashed line), and supercritical states.Fig. 4. Stochastic resonance.a Phase diagram of stochastic resonance. The color represents the normalized time-averaged phase lag derivative, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }\rangle /{\omega }_{\eta }$$\end{document} . The vertical axis, ωη/ωa, quantifies the relative intensity of turbulent fluctuations (noise) to the magnetic field strength. The horizontal axis, ωH/ωa, distinguishes subcritical, critical (ωH/ωa = 1/2, marked by the black dot-dashed line) and supercritical regimes. For a fixed ωH/ωa in the subcritical (supercritical) range, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }\rangle /{\omega }_{\eta }$$\end{document} exhibits a non-monotonic (monotonic) dependence on noise. A pronounced resonant peak appears at ωη/ωa ≈ 1 when ωH/ωa < 1/2 (b), with a perfect collapse in the linear regime when normalized by ωH/ωa (b, inset). The diamond symbols in (a) denote the parameter settings for experiments shown in (c). c Experimental evidence of stochastic resonance (diamonds) at ωH/ωa = 0.12 ± 0.01. For comparison, simulation results at a similar value of ωH/ωa (circles) are also shown. The error bars correspond to the standard deviation calculated from 20 independent measurements at the same parameter value. d Effective potential well framework: the magnetic field contributes a cosine potential modulated by a linear term, and turbulence acts as a stochastic noise, see Eq. (2). When ωH/ωa < 1/2 (e), for weak noise, the phase lag remains near the local minimum of the potential well. As noise increases, the particle phase lag exhibits unidirectional escape, signaling stochastic resonance. At even higher noise levels, escape becomes bidirectional, leading to a subsequent decline in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }\rangle$$\end{document} . For ωH/ωa ≪ 1/2, the potential well is sufficiently deep to strongly localize the particle phase lag (f). At the critical threshold ωH/ωa = 1/2, the potential well vanishes, making the phase lag sensitive to small perturbations (g). In the supercritical range (ωH/ωa ≫ 1/2), the potential well inherently favors persistent phase slipping, regardless of noise intensity (h).
For a fixed ωH/ωa in the subcritical regime (ωH/ωa < 1/2), the particle rotation in the co-rotating frame \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }\rangle /{\omega }_{\eta }$$\end{document} exhibits a non-monotonic dependence on noise intensity (Fig. 4b), with a pronounced resonant peak at ωη/ωa ≈ 1, indicating that resonance arises at the transition from the “phase-locked” regime to the “turbulent-dynamics” regime. This is a clear suggestion that the observed resonance behavior results from the cooperation between turbulent fluctuations and magnetic forcing. SR is also clearly observed experimentally, as shown in Fig. 4c: diamonds denote the experimental results at ωH/ωa = 0.12 ± 0.01 (corresponding to the points marked in the phase diagram in Fig. 4a), while circles represent the simulation results at a comparable value of ωH/ωa for direct comparison. The experimental and simulation results agree reasonably well; at the resonance peak, the experimental curve appears flatter which can be attributed to effects of polydispersity of the particles used in experiments. The resonance mechanism can be interpreted via an effective potential well framework (Fig. 4d), where the magnetic field induces a cosine-type potential (i.e., \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$-\frac{{\omega }_{{{{\rm{a}}}}}}{4}\cos 2\beta$$\end{document} ) modulated by a linear term (i.e., −ωHβ), while turbulence acts as a stochastic noise.
For ωH/ωa < 1/2 (Fig. 4e), the phase lag angle is initially trapped in the potential minimum at weak noise levels (small ωη/ωa), exhibiting only mild fluctuations locally. As the turbulence-induced fluctuations become comparable to the magnetic field strength (ωη/ωa ≈ 1), preferential unidirectional escape occurs. This clearly identifies the observed resonant behavior as an example of SR: a cooperative interplay between stochastic forcing and the external deterministic forcing. This phenomenon can be thought of as the particle “surfing” the rotating magnetic field, aided by the turbulent vorticity. At even higher noise levels, randomization overrides the cooperative effect, causing bidirectional escape and a decline in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }\rangle$$\end{document} .
At the critical threshold ωH/ωa = 1/2, the potential well vanishes, rendering β critically sensitive to any perturbations (Fig. 4g). In the supercritical range (ωH/ωa ≫ 1/2), the potential well inherently favors persistent phase slipping, regardless of noise intensity (Fig. 4h), such that the SR mechanism vanishes in this regime.
In strongly subcritical regime (ωH/ωa ≪ 1/2), the resonance peak is purely set by the turbulent noise, making it sensitive to the preferential sampling effect^7,9,10^ (see Supplementary Information for neutral and heavy particles results). Here we introduce the relaxation time \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\tau }^{*}=1/\left({\omega }_{{{{\rm{a}}}}}\sqrt{1-4{\left(\frac{{\omega }_{H}}{{\omega }_{{{{\rm{a}}}}}}\right)}^{2}}\right)$$\end{document} . Physically, τ^^ represents the characteristic time it takes for a perturbed particle rotational motion to return to its phase-locked state, i.e., \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{\beta }=0$$\end{document} , after a small disturbance (detailed derivations are in Supplementary Information). Close to the critical ratio ω_H_/ωa = 1/2, the significantly increasing particle relaxation time τ^^ shifts the resonance position, making it the dominant factor. Therefore, a longer τ^*^ allows particles to escape the potential well even at lower noise levels.
Symmetry breaking mechanism
A striking consequence of SR is the symmetry-breaking effect, manifesting as a “zero plus zero is greater than zero” phenomenon. This leverages the fact that in the weak noise regime (ωη < ωa), the response of the system to the noise is non-linear (Fig. 4b, purple-shaded part).
To reveal this effect, instead of applying a constant-frequency rotating magnetic field (ωH), we now impose a time-dependent field rotation with a zero-average, ωH(t) (Fig. 5a). Each oscillation cycle consists of a forward and a reverse rotation phase, characterized by angular velocities \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\omega }_{H}^{+}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\omega }_{H}^{-}$$\end{document} with corresponding durations \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${T}_{H}^{+}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${T}_{H}^{-}$$\end{document} . These parameters are scaled as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\omega }_{H}^{+}/{\omega }_{H}^{-}={T}_{H}^{-}/{T}_{H}^{+}=\gamma$$\end{document} , ensuring that the net angular velocity over a full cycle remains zero. Two distinct response regimes emerge as ωη/ωa increases for a fixed γ (Fig. 5b), highlighting the ability of this oscillation framework to probe the nonlinearity of the system in response to noise. Specifically, the response \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }\rangle /{\omega }_{\eta }$$\end{document} decreases monotonically in the nonlinear regime (weak noise limit), while it remains nearly constant and close to zero in the linear regime (strong noise limit, characterized by response collapse, see also Fig. 4b, inset). This behavior aligns well with a superposition model, which predicts the response as
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }\rangle=\frac{\langle \dot{\beta }({\omega }_{H}^{+})\rangle {T}_{H}^{+}-\langle \dot{\beta }({\omega }_{H}^{-})\rangle {T}_{H}^{-}}{{T}_{H}},$$\end{document}with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${T}_{H}={T}_{H}^{+}+{T}_{H}^{-}$$\end{document} , where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }({\omega }_{H}^{+})\rangle$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }({\omega }_{H}^{-})\rangle$$\end{document} are interpolated from the numerical results in Fig. 4b. The predicted trends (Fig. 5b, dashed lines) show excellent agreement with the frequency-modulated field results (Fig. 5b, symbols).Fig. 5. Symmetry breaking.a The applied magnetic field oscillates with zero-average angular velocity, ωH(t). Without loss of generality, we set \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\omega }_{H}^{+}/{\omega }_{H}^{-}={T}_{H}^{-}/{T}_{H}^{+}=\gamma \ge 1$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\omega }_{H}^{+}=0.9{\omega }_{{{{\rm{cr}}}}}$$\end{document} . b The normalized phase lag derivative, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }\rangle /{\omega }_{\eta }$$\end{document} , shows distinct regimes as noise intensity ωη/ωa increases for fixed γ. In the nonlinear regime (low noise), \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\langle \dot{\beta }\rangle /{\omega }_{\eta }$$\end{document} decreases monotonically, while in the linear regime (strong noise, where system responses collapse, see Fig. 4b, inset), it remains nearly unchanged and close to zero. When γ = 1, the particle rotation motions cancel each other out in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${T}_{H}^{-}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${T}_{H}^{+}$$\end{document} , resulting in a zero response. The observed deviation from a strictly zero value in the simulations is attributable to the finite periods of the applied oscillation cycles of ωH(t). Theoretically, with an infinite simulation time, the value would converge to zero. The error bars represent the standard deviation of the results obtained by sampling over one full period of ωH(t) at the final stage of the simulation. The predicted trends by the superposition model of Eq. (3) (b, dashed lines) exhibit excellent agreement with the numerical data (b, symbols). c Particle motion in laboratory frame. With turbulence alone (green dotted line), the particle shows no net angular movement in terms of ψ; with only the magnetic field (blue dashed line), the motion is periodic with zero net angular displacement. When both effects act together, a net angular movement emerges (black solid line), revealing symmetry-breaking mechanism. This is different from the case of constant ωH (inset), where the particle follows ψ = ωHt (black solid line) when it is purely driven by a magnetic field. Parameter settings are TH = 50τη ≈ TL with TL the integral time scale of turbulence, ωη/ωa = 0.25, and γ = 2. Inset: ωη/ωa = 0.25 and ωH/ωa = 0.4.
We emphasize that when combining turbulence and the frequency-modulated magnetic field in this way, this leads to particle motion directly in the lab frame. To evidence this, we consider the lab-frame angular particle motion ψ(t) (Fig. 5c). With turbulence alone, the particle exhibits no net angular movement (green dotted line). Under a purely frequency-oscillatory magnetic field, the particle undergoes periodic motion with zero net angular displacement (blue dashed line). However, when both effects act concurrently, a net, non-zero mean angular drift remarkably emerges (black solid line), as the turbulent vorticity aids the particle in “surfing” the magnetic field, which reveals an unexpected symmetry-breaking mechanism in turbulent flows. This stands in contrast to the case of a constant ωH (Fig. 5c, inset), where a purely magnetic field-driven particle, instead of showing zero net angular movement as in the oscillating ωH scheme, follows ψ = ωHt (black solid line).
Discussion
We have revealed the first experimental and numerical evidence of SR for magnetic particles suspended in turbulence, demonstrating how the interplay between turbulent fluctuations and external forcing induces coherent rotational motion of magnetic particles. Three distinct regimes: “phase-locked”, “back-and-forth”, and “turbulent-dynamics” are identified, with SR occurring at the transition between “phase-locked” and “turbulent-dynamics” regimes. Remarkably, even when both turbulence and the magnetic field have zero mean angular velocity, their collective effect is capable of generating a nonzero particle rotational response, as the turbulent vorticity helps the particle to “surf” the rotating magnetic field.
Different particle properties influence their response, and when combined with the SR mechanism, they enable a novel magnetic resonance-based approach to measure turbulent vorticity. This method leverages activated magnetic particles as probes, effectively functioning as a turbulent vorticity “microscope” that operates in various conditions. In the present work, we validate this concept using optical measurements of particle rotation. However, in principle, once the angular response is well-characterized, the same protocol could be applied even in optically inaccessible environments by measuring the magnetic fields emitted by the spinning magnetic particle directly. This concept shares similar principles with Magnetic Particle Spectroscopy, a technique that uses the dynamic magnetic responses of magnetic particles to remotely sense properties of their surrounding environment^54^. In such scenarios, the activated magnetic particles would function as remotely readable vorticity probes. The protocol is straightforward: by dispersing magnetic particles into turbulence and varying the magnetic frequency ωH and intensity ωa while maintaining a constant ratio ωH/ωa, one can measure the resonance curve of the averaged particle angular response. The location of the resonance peak directly reveals the turbulent vorticity magnitude.
While the present work focuses on understanding particle angular dynamics, it lays essential groundwork for subsequent investigations into how these actively tunable particles can shape turbulent structures, opening possibilities for flexible manipulation of turbulence, where propeller-like particles could be tuned to navigate specific flow structures, inject energy into turbulence, and influence flow properties such as intermittency^55^. More broadly, our study provides fundamental mechanisms that externally tunable rotational dynamics could be leveraged to design active materials^56–58^. If organized coherently, their collective behavior may also exhibit hydrodynamic properties analogous to systems with odd viscosity, which may be promising experimental realizations of chiral fluids^59,60^.
This research paves the way for innovative experimental techniques, offering a simple yet underutilized approach to studying and controlling turbulence dynamics.
Methods
Experimental setup
The experimental setup comprises three fundamental components: a turbulence generator, a magnetic field generator, and magnetic particles. Detailed specifications of the experimental apparatus, along with a comprehensive description of the particle rotation tracking methodology, are not relevant to the core physics of current work, and hence they are provided in our dedicated instrumentation and methodology paper^61^. Here, we present only the essential information pertinent to the current study.
Turbulence generator: The experiments are conducted in a confined Von Kármán-type turbulent flow, which serves as a model system for studying HIT, particularly in the context of Eulerian and Lagrangian statistics^48,49,62–64^.
The working fluid is water, contained within an octagonal-shaped vessel with an internal diameter of 2R = 150 mm and a height of 220 mm. This specific geometry is chosen to enhance optical accessibility for visualization. The flow is driven by two counter-rotating disks, each with a diameter of 150 mm, separated by a distance of 150 mm, and rotating at the same angular velocity magnitude in opposite directions. Each disk is equipped with six straight blades, each 10 mm in height, designed to enhance inertial stirring. The disks are powered by two calibrated DC motors operating at a constant voltage, ensuring that the angular velocity Ω remains stable over time with a precision of approximately 1%. The motor rotation frequency can be adjusted within the range of 0.25–8 Hz. Under these conditions, a fully developed turbulent flow with up to Re_λ_ = 447 can be generated in the central region of the container with small-scale statistics close to local isotropy^46–49,62–64^. This region, measuring approximately 20 × 20 × 20 mm, serves as the primary domain for measurements. To enable particle rotational motion tracking, a Photron NOVA S6 high-speed camera, operating at 3000 frames per second with a resolution of 1024 × 1024 pixels, is employed to record particle motion within the flow field.
Magnetic field generator: A uniform rotating planar magnetic field within a restricted measurement volume, H = Hh, with an angular frequency ωH, is generated using a system of two pairs of mutually perpendicular Helmholtz coils. The vertical pair (marked in red in Fig. 1a) has an interior diameter of 300 mm, while the horizontal pair (marked in blue in Fig. 1a) has an interior diameter of 245 mm. For each pair, the separation between the two coils is equal to the corresponding coil diameter.
Each Helmholtz coil pair generates a magnetic field along a single axis. By driving the two perpendicular coil pairs with AC currents that have a phase difference of π/2, a resultant magnetic field is produced that rotates at a constant angular velocity ωH. This system is capable of generating rotating magnetic fields with frequencies of up to 50 Hz.
Magnetic particles: The magnetic particles are produced by coating spherical Styrofoam cores with a layer of magnetic paint, which is paramagnetic. The Styrofoam cores are initially arranged in a fixed, evenly spaced configuration on a substrate. The magnetic paint is then applied from the top of the substrate to the Styrofoam spheres through a spraying process, ensuring that all particles are coated in a consistent manner; however, at the single-particle level the coating is not perfectly homogeneous, and small variations in thickness or distribution lead to slight differences in magnetic properties between particles. Finally, the coated particles are detached from the substrate.
The particle volume fraction used in the experiments is 0.17%, ensuring that the system remains in the dilute regime. This low concentration minimizes inter-particle interactions, allowing us to focus on the dynamics of particle response to the complex interplay between turbulent fluctuations and the external magnetic field, without significant effects from collective dynamics.
These particles exhibit anisotropy in three key aspects: (1) Shape: The particles are approximately spherical, with a mean diameter of 0.762 ± 0.066 mm. (2) Density: The magnetic particles have a mean density of ρparticle = 208 ± 14 kg/m^3^. (3) Magnetic Properties: Due to the non-uniform distribution of the magnetic paint on the particle surface, the magnetic properties vary between particles. Specifically, the magnetic moment magnitude ∣m∣, the orientation of the magnetic moment, and the spatial distribution of magnetization across the surface differ from particle to particle. This variability influences how individual particles respond to an external magnetic field. All these three aspects of particle polydispersity can introduce minor differences between the simulation results and the experimental observations. However, these differences remain small, and overall, the simulations show good agreement with the experiments. More importantly, our theoretical modeling and the simulations successfully capture key physical features observed in the experiments, such as the location of peak responses.
Particle rotation tracking: The particle rotational motion is tracked using an efficient particle rotation tracking algorithm. The algorithm initiates by pre-processing the experimental images through the application of Gaussian filtering to mitigate noise, contrast enhancement to optimize image quality, and central cropping to ensure precise target particle localization. Subsequently, a set of discrete candidate rotation angles is defined, and corresponding theoretical rotated images are generated utilizing 3D rotation matrices. For each temporal frame, the algorithm computes the cross-correlation coefficient between the actual experimental image and each theoretical rotated image to determine the optimal matching rotation angle. By tracking these optimal angles across the time series, the rotational motion trajectory of the particle is reconstructed and visualized. More details of the tracking methodology are reported in our dedicated instrumentation and methodology paper^61^.
Theoretical model
A paramagnetic particle moves in a homogeneous and isotropic turbulent flow while subjected to an external magnetic field H. The magnetic field rotates in the xoy-plane (laboratory frame) with an angular velocity ωH and intensity H, expressed as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{\bf{H}}}}=H{{{\bf{h}}}}=H{(\cos ({\omega }_{H}t),\sin ({\omega }_{H}t),0)}^{{{{\rm{T}}}}}$$\end{document} . The particle is modeled as a sphere with intrinsic magnetic anisotropy (same as that of a prolate spheroid of aspect ratio 1.2), characterized by a preferred magnetization direction n, which dynamically evolves over time. This instantaneous unit vector is given by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{\bf{n}}}}={(\sin \alpha \cos \psi,\sin \alpha \sin \psi,\cos \alpha )}^{{{{\rm{T}}}}}$$\end{document} , where α and ψ (defined in Fig. 2a) describe the orientation of n.
When the particle is placed in an external magnetic field, H, its induced magnetic moment takes the form^19^
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{\bf{m}}}}={{{{\mathcal{V}}}}}_{{{{\rm{p}}}}}H[\,{\chi }_{\perp }{{{\bf{h}}}}+({\chi }_{\parallel }-{\chi }_{\perp })({{{\bf{n}}}}\cdot {{{\bf{h}}}}){{{\bf{n}}}}],$$\end{document}where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\mathcal{V}}}}}_{{{{\rm{p}}}}}=4\pi {r}^{3}/3$$\end{document} is the particle volume, with r denoting its radius. The anisotropic magnetic susceptibilities χ∥,⊥, arising due to the anisotropy of the particle shape (e.g., spheroidal) and the demagnetizing field factors N∥,⊥, and χ∥,⊥ = χ0/(1 + χ0N∥,⊥)^51^, causes the particle to respond differently to an external magnetic field depending on the direction. The subscripts parallel (∥) and perpendicular (⊥) refer to directions relative to the preferred magnetization direction n. Due to the tendency of the magnetic moment m to align with the magnetic field direction h, the rotating magnetic field induces the particle rotational motion. The magnetic torque acting on the particle is given by TM = μwm × H, where μw is the magnetic permeability. The magnetic torque tries to align the particle’s preferred magnetic axis (n) with the magnetic field direction (h) and is balanced by the rotational drag torque (approximated in the Stokes regime) TD = −ξr(ωp–ωf)^65^, where ξr = 8 πμr^3^ is the rotational drag coefficient and μ is the dynamic viscosity of water.
The particle rotational velocity ωp naturally decomposes into two components: a perpendicular component ωp,⊥ (tumbling) and a parallel component ωp,∥ (spinning) along n.
The particle rotational motion is governed by the balance of torques, and we consider the overdamped condition^50^
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\bf{T}}}}}_{{{{\rm{M}}}}}+{{{{\bf{T}}}}}_{{{{\rm{D}}}}}={{{\mathbf{0}}}},$$\end{document}with TM the magnetic torque and TD the rotational drag torque.
By projecting this balance along and perpendicular to n, we obtain the equations governing tumbling and spinning. The tumbling dynamics, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{{{{\bf{n}}}}}$$\end{document} , obey
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\frac{d{{{\bf{n}}}}}{dt}={\omega }_{{{{\rm{a}}}}}({{{\bf{n}}}}\cdot {{{\bf{h}}}})\left({{{\bf{h}}}}-({{{\bf{n}}}}\cdot {{{\bf{h}}}}){{{\bf{n}}}}\right)+{{{{\boldsymbol{\omega }}}}}_{{{{\rm{f}}}}}\times {{{\bf{n}}}},$$\end{document}where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\omega }_{{{{\rm{a}}}}}=\mu {{{{\mathcal{V}}}}}_{{{{\rm{p}}}}}{H}^{2}({\chi }_{\parallel }-{\chi }_{\perp })/{\xi }_{{{{\rm{r}}}}}$$\end{document} is the normalized magnetic field intensity, which is also called anisotropic frequency^19,50^, with ωa > 0 (ωa < 0) for a prolate (oblate) particle.
The particle spinning rate ωp,∥ satisfies
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\frac{d\left({\xi }_{{{{\rm{r}}}}}{{{{\boldsymbol{\omega }}}}}_{{{\rm{p}}},\parallel }\right)}{dt}=-{\xi }_{{{{\rm{r}}}}}\left({{{{\boldsymbol{\omega }}}}}_{{{\rm{p}}},\parallel }-\left({{{{\boldsymbol{\omega }}}}}_{{{{\rm{f}}}}}\cdot {{{\bf{n}}}}\right){{{\bf{n}}}}\right),$$\end{document}where for the overdamped particle, we have \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\frac{d\left({\xi }_{{{{\rm{r}}}}}{{{{\boldsymbol{\omega }}}}}_{{{\rm{p}}},\parallel }\right)}{dt}={{{\mathbf{0}}}}$$\end{document} . ωf is half of the vorticity of the background flow field, which is obtained through performing Direct Numerical Simulations (DNS) using a pseudo-spectral method. The validity of the overdamped model is further discussed in the Supplementary Information.
The total particle angular velocity is thus \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\boldsymbol{\omega }}}}}_{{{{\rm{p}}}}}={{{{\boldsymbol{\omega }}}}}_{{{\rm{p}}},\parallel }+{{{{\boldsymbol{\omega }}}}}_{{{\rm{p}}},\perp }=({{{\bf{n}}}}\cdot {{{{\boldsymbol{\omega }}}}}_{{{{\rm{p}}}}}){{{\bf{n}}}}+{{{\bf{n}}}}\times \dot{{{{\bf{n}}}}}$$\end{document} , with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{{{{\bf{n}}}}}=d{{{\bf{n}}}}/dt$$\end{document} (see Eq. (6)). Equations of motion (6) and (7) are numerically integrated using a third-order Runge–Kutta method.
The governing Eq. (5) can be reformulated in terms of the phase lag angle β and the orientation angle α, where α denotes the angle between the unit vector n and the z-axis, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hat{{{{\bf{z}}}}}$$\end{document} .
To derive the governing equation for β, we project Eq. (5) onto \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hat{{{{\bf{z}}}}}$$\end{document} -direction, yielding
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\bf{T}}}}}_{{{{\rm{M}}}}}\cdot \hat{{{{\bf{z}}}}}+{{{{\bf{T}}}}}_{{{{\rm{D}}}}}\cdot \hat{{{{\bf{z}}}}}=0.$$\end{document}Here, the magnetic torque component along \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hat{{{{\bf{z}}}}}$$\end{document} is given by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\bf{T}}}}}_{{{{\rm{M}}}}}\cdot \hat{{{{\bf{z}}}}}={\omega }_{{{{\rm{a}}}}}{\xi }_{{{{\rm{r}}}}}{\sin }^{2}\alpha \sin \beta \cos \beta$$\end{document} , while the drag torque contribution is expressed as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\bf{T}}}}}_{{{{\rm{D}}}}}\cdot \hat{{{{\bf{z}}}}}=-{\xi }_{{{{\rm{r}}}}}({\omega }_{{{\rm{p}}},z}-{\omega }_{{{\rm{f}}},z})$$\end{document} .
Rewriting Eq. (8) explicitly in terms of β and α, we obtain
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{\beta }={\omega }_{H}-\frac{{\omega }_{{{{\rm{a}}}}}}{2}\kappa \sin 2\beta -{\omega }_{{{\rm{f}}},z},$$\end{document}with the coefficient \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\kappa={\sin }^{2}\alpha$$\end{document} .
To establish the evolution equation for α, we introduce a unit vector \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{\bf{N}}}}={(-\sin \psi,\cos \psi,0)}^{{{{\rm{T}}}}}$$\end{document} , which denotes the direction perpendicular to nxoy, the projection of n onto the xoy-plane (see Fig. 2a for the definition). Projecting Eq. (5) along N leads to
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\bf{T}}}}}_{{{{\rm{M}}}}}\cdot {{{\bf{N}}}}+{{{{\bf{T}}}}}_{{{{\rm{D}}}}}\cdot {{{\bf{N}}}}=0.$$\end{document}Here, the magnetic torque component along N is given by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\bf{T}}}}}_{{{{\rm{M}}}}}\cdot {{{\bf{N}}}}=\frac{1}{2}{\omega }_{{{{\rm{a}}}}}{\xi }_{{{{\rm{r}}}}}{\cos }^{2}\beta \sin 2\alpha$$\end{document} , while the drag torque contribution is TD ⋅ N = −ξr(ωp ⋅ N − ωf ⋅ N). Noting that \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\boldsymbol{\omega }}}}}_{{{{\rm{p}}}}}\cdot {{{\bf{N}}}}=\dot{\alpha }$$\end{document} , Eq. (10) can be rewritten explicitly as
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{\alpha }={{{{\boldsymbol{\omega }}}}}_{{{{\rm{f}}}}}\cdot {{{\bf{N}}}}+\frac{{\omega }_{{{{\rm{a}}}}}}{2}\sin 2\alpha \,\,{\cos }^{2}\beta .$$\end{document}To capture the essential physics of the system, we assume κ ≈ 1, which is reasonable because κ follows a probability density function that is sharply peaked at 1 and decays rapidly to zero for values smaller than unity within the investigated parameter space (further details are provided in the Supplementary Information). The equation for β can be reduced to
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{\beta }={\omega }_{H}-\frac{{\omega }_{{{{\rm{a}}}}}}{2}\sin 2\beta -{\omega }_{{{\rm{f}}},z},$$\end{document}which is controlled closely by the magnetic field.
Equation (12) is used to predict the phase transition between the phase-locked regime and the turbulent-dynamics and back-and-forth regimes (blue solid line in Fig. 3m), and is ultimately describing the dynamics that gives rise to the mechanism of SR (Fig. 4c–g). Similar dynamical equations as Eq. (12) have been widely applied in multiple physical systems, including condensed matter physics, chemical and biological systems, optics and laser physics, electrical engineering, and financial econometrics^66,67^, with a range of related physical problems, such as Josephson junctions^68–72^, vortex diffusion and flux-line dynamics in superconductors^68,73,74^, reaction-rate characteristics relevant to chemistry, engineering, and biology^75^, Brownian motors^76^, and decision-making processes under uncertainty which resemble financial market fluctuations^67,72^. The universality of this equation suggests that our findings have broad implications beyond the specific context of magnetic particles in turbulence, offering deeper theoretical insights into nonlinear dynamical systems and facilitating the development of new applications in diverse scientific and technological domains.
Now we show the details of predicting the boundary (Fig. 3m, blue solid line) that separates the phase-locked regime (Fig. 3m, purple region) from the back-and-forth and turbulent-dynamics regimes. The phase-locked regime exists if the equation \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{\beta }=0$$\end{document} admits fixed points. Based on the simplified governing equation, Eq. (2) or Eq. (12), this condition translates to ensuring that the following equation has real solutions for any imposed ωH:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\frac{{\omega }_{{{{\rm{a}}}}}}{2}\sin 2\beta={\omega }_{H}-{\omega }_{{{{\rm{f}}}},z}.$$\end{document}The left-hand side of Eq. (13) is bounded within [−ωa/2, ωa/2], while the right-hand side varies typically over \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\left[{\omega }_{H}-\sqrt{\sigma ({\omega }_{{{\rm{f}}},z})},{\omega }_{H}+\sqrt{\sigma ({\omega }_{{{\rm{f}}},z})}\right]$$\end{document} , where σ(ωf,z) is the variance of ωf,z. For Eq. (13) to have real solutions for any imposed ωH, the following condition must hold:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\omega }_{H}+\sqrt{\sigma ({\omega }_{{{\rm{f}}},z})}\le \frac{{\omega }_{{{{\rm{a}}}}}}{2},$$\end{document}which gives
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\omega }_{{{{\rm{a}}}}}\ge 2{\omega }_{H}+2\sqrt{\sigma ({\omega }_{{{\rm{f}}},z})}.$$\end{document}Notice that the standard deviation of ωf,z is approximated to be of the order of the Kolmogorov frequency scale, i.e., \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sqrt{\sigma ({\omega }_{{{\rm{f}}},z})} \sim {\omega }_{\eta }$$\end{document} , where ωη ~ 1/τη. In our simulations, we observe that \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$2\sqrt{\sigma ({\omega }_{{{\rm{f}}},z})} \sim {\omega }_{\eta }$$\end{document} , which justifies this approximation. Based on this, we obtain the final criterion for the phase-locked regime from Eq. (15) as
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\omega }_{{{{\rm{a}}}}}\ge 2{\omega }_{H}+{\omega }_{\eta }.$$\end{document}Direct numerical simulation of the Navier-Stokes equation
To accurately integrate the governing Eqs. (6) and (7), it is essential to correctly incorporate the Lagrangian vorticity signal of the particles ωf.
We begin by performing DNS of statistically stationary homogeneous isotropic turbulence (HIT) to obtain the Lagrangian vorticity dynamics (ωf) of small, passively advected particles. The computational domain consists of a cubic box of size L = 2π with periodic boundary conditions. The simulations employ a pseudo-spectral method^52^ for spatial discretization and an Adams-Bashforth scheme for temporal integration, with the standard 2/3 dealiasing rule applied to ensure numerical accuracy^53^. The numerical framework has been validated against multiple integration schemes, interpolation methods, and large-scale forcing techniques^77^.
The suspended particles are modeled as point-like, dilute tracers that interact with the turbulent flow through hydrodynamic forces alone. Their dynamics are characterized by the Stokes number (St) and the density contrast (βp), which describe their response time relative to the fluid and their mass relative to the surrounding medium, respectively. To capture the essential physics while simplifying the problem, we adopt a one-way coupled model, which has been extensively validated in both experimental and numerical studies^9,10,78–82^.
The governing equation for particle motion is
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\ddot{{{{\boldsymbol{x}}}}}}_{{{{\rm{p}}}}}={\beta }_{{{{\rm{p}}}}}\frac{{{{\rm{D}}}}{[{{{{\bf{u}}}}}_{{{{\rm{f}}}}}]}_{{{{\rm{O}}}}}}{{{{\rm{D}}}}t}-\frac{1}{\,{{\rm{St}}}\,}\left({\dot{{{{\bf{x}}}}}}_{{{{\rm{p}}}}}-{[{{{{\bf{u}}}}}_{{{{\rm{f}}}}}]}_{{{{\rm{O}}}}}\right),$$\end{document}where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${[{{{{\bf{u}}}}}_{{{{\rm{f}}}}}]}_{{{{\rm{O}}}}}$$\end{document} denotes the fluid velocity with the subscript O indicating evaluation at the center of the particle, xp represents the particle position, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\bf{v}}}}}_{{{{\rm{p}}}}}={\dot{{{{\bf{x}}}}}}_{{{{\rm{p}}}}}$$\end{document} is the particle velocity. The Lagrangian vorticity is computed as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{{{\boldsymbol{\omega }}}}}_{{{{\rm{f}}}}}=\frac{1}{2}{[{{{\boldsymbol{\nabla }}}}\times {{{{\bf{u}}}}}_{{{{\rm{f}}}}}]}_{{{{\rm{O}}}}}$$\end{document} . The density contrast is defined as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\beta }_{{{{\rm{p}}}}}=\frac{3}{1+2{\rho }_{{{{\rm{p}}}}}/{\rho }_{{{{\rm{f}}}}}}$$\end{document} , where ρp/ρf is the particle-to-fluid density ratio. The Stokes number is given by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\,{{\rm{St}}}={d}_{{{{\rm{p}}}}}^{2}/(12{\beta }_{{{{\rm{p}}}}}\nu {\tau }_{\eta })$$\end{document} , where dp is the particle diameter and ν is the fluid kinematic viscosity.
In this study, we perform three-dimensional DNS at a resolution of N^3^ = 128^3^ with a Taylor-scale Reynolds number Re_λ_ ≈ 60^10,52,53,77^ and collect 512 independent samples for three particle types: light (βp → 3.0, representing bubbles where ρp → 0), neutral (βp = 1.0), and heavy (βp = 0.01) particles, all with a Stokes number of St = 1.0. This specific Stokes number is chosen because light particles in HIT tend to accumulate preferentially in vortex filaments, leading to strong dynamic localization^7,9,10,83^.
A key consideration in extending our results to higher Reynolds number turbulence is the effect of intermittency, particularly the heavy-tailed distribution of vorticity fluctuations in fully developed turbulence. In principle, for very large Reynolds numbers, the intermittency of the small-scale vorticity field should be accounted for, as extreme vorticity events become more pronounced. However, previous studies have shown that while intermittency affects the fine-scale structure of turbulence, its influence on the rotational dynamics of small particles (compared to the Kolmogorov length scale), particularly in the parameter regimes relevant to our study, remains relatively weak, because the intense vorticity events, although frequent, are spatially localized and short-lived; thus, their cumulative effect on particle rotation is averaged out^81,84,85^. Our results and the underlying mechanisms in the main paper provide valuable insights into the fundamental principles governing particle rotation in turbulent flows. Future research could explore the effects of finite-sized particles and beyond the dilute regime.
Supplementary information
Supplementary information Description of additional Supplementary Files Supplementary Movie 1 Supplementary Movie 2 Supplementary Movie 3 Supplementary Movie 4 Supplementary Movie 5 Supplementary Movie 6 Supplementary Movie 7 Supplementary Movie 8 Supplementary Movie 9 Supplementary Movie 10 Transparent Peer Review file
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1Frisch, U. Turbulence: The Legacy of AN Kolmogorov (Cambridge University Press, 1995).
- 2Pope, S. B. Turbulent Flows (Cambridge University Press, 2000).
- 3Guazzelli, E. & Morris, J. F. A Physical Introduction to Suspension Dynamics, Vol. 45 (Cambridge University Press, 2011).
- 4Peyret, R. Spectral Methods for Incompressible Viscous Flow, Vol. 148 (Springer, 2002).
- 5Freitas, A., de Wit, X. M., Wang, Z., Biferale, L. & Toschi, F. Statistical properties of turbulence under a smart Lagrangian forcing. Preprint at https://arxiv.org/abs/2508.06660 (2025).
- 6Wu, C. et al. Tracking the rotation of light magnetic particles in turbulence. Preprint at https://arxiv.org/abs/2506.21769 (2025).
- 7Coffey, W. & Kalmykov, Y. P. The Langevin Equation: With Applications to Stochastic Problems in Physics, Chemistry and Electrical Engineering, Vol. 27 (World Scientific, 2012).
- 8Calzavarini, E. et al. Effect of Faxén forces on acceleration statistics of material particles in turbulent flow. In Advances in Turbulence XII (ed. Eckhardt, B.) 11–14 (Springer Berlin Heidelberg, 2009).
