## INTRODUCTION

The propagation of elastic waves in soft materials plays a crucial role in the spatiotemporal transmission of mechanical signals, e.g., in biological mechanotransduction (*1*, *2*) or in the failure of marginal solids (*3*–*6*). At high Reynolds numbers *Re* ≫ 1, inertia dominates and wave propagation can be readily observed (*7*–*9*). However, mechanical cues in soft and biological materials often occur at low *Re* (*10*), playing a key role in their spatiotemporal response to mechanical perturbations (*11*, *12*). For example, in marginally stable systems, such as jammed packings or fiber networks, a stress at the right position can cause total loss of rigidity (*3*–*6*). In living cells, the propagation of mechanical signals through soft structures is crucial in mechanotransduction (*1*, *2*) and controls, for example, cell differentiation (*13*). Also, in soft robotics, sensing and control require the transmission of mechanical signals over a distance (*14*).

Localized mechanical signals can spread in space and time by propagation of elastic waves (*7*–*9*). This wave propagation is intimately linked to the mechanical properties of the medium (*15*). In very soft materials, the transmission of mechanical signals poses challenges, because their intrinsically dissipative nature leads to energy losses as the wave propagates (*10*, *16*). In particular, at low excitation frequencies, corresponding to low Reynolds numbers, viscous attenuation of the wave signal is strong, and its detection is challenging. Moreover, in ultrasoft solids, the relative amplitudes of thermal fluctuations are large, thus further obscuring robust wave propagation and their detection in experiments.

In this Letter, we show how Fourier filtering can reveal even very weak and strongly damped elastic waves at extremely low Reynolds numbers, *Re* ∼ 10^{−6}, in ultrasoft solids formed from crystals and glasses of colloids in two dimensions. We create a localized oscillatory perturbation within these solids with an optical tweezer and use video microscopy and frequency domain filtering to quantify the spatiotemporal strain response. On the basis of an overdamped equation of motion, which is in excellent quantitative agreement with our experimental results, this enables a full characterization of low *Re* wave propagation and attenuation in ultrasoft solids. Moreover, we show how the analysis of these waves can be used to obtain the linear mechanical moduli of very weak elastic solids.

## RESULTS

We prepare two-dimensional (2D) colloidal crystals by sedimenting monodisperse silica particles with radius *a* = 3.04 μm suspended in an aqueous buffer at pH 8.4 and an ionic strength of 10 mM, giving a Debye screening length κ^{−1} ≈ 3 nm ≪*a*. This short-ranged electrostatic repulsion yields dense 2D crystals with long-ranged hexagonal order at a packing fraction of 0.84 (Fig. 1A and section S1.1).

To create a propagating mechanical wave, we trap a single particle of the crystal in an optical trap and force it into an in-plane oscillatory motion with an amplitude *A*_{0} = 2.5 μm (Fig. 2A). We vary the frequency ω of this motion between 0.05 and 10 rad/s, corresponding to Reynolds numbers between 4 · 10^{−7} and 8 · 10^{−5} (*17*), and confirm that this mechanical excitation is well within the linear regime (section S2).

Note that the self-diffusion time of particles in the crystal, τ = *a*^{2}/4*D* ≈ 100 s, so that the Deborah number *De* = ωτ is larger than 1 in all our experiments (section S3), implying that the material can be considered a viscoelastic solid.

The oscillating particle creates a mechanical wave that propagates through the surrounding material. This sets up a net ballistic displacement of the particles adjacent to the oscillating bead. However, the signal of interest is convolved with the inherent Brownian motion of these microscopic colloids. In particular, far away from the trapped colloid, where the elastic signal is attenuated, it may drown in the Brownian noise (Fig. 2, B and E). Upon filtering the positional trajectory of each particle in the frequency domain at the driving frequency (Fig. 2C), even small displacements due to the propagating wave become apparent (Fig. 2D). To ensure statistical reliability of these data, we set up a real-time distributed particle tracking algorithm that allows us to collect data during 20,000 to 35,000 frames, which is equivalent to 50 to 500 oscillation cycles. While the unfiltered mean square displacement of the particles shows no apparent signature of the perturbation, except at the forced particle (Fig. 2E), the Fourier-filtered amplitude map shows a propagating mechanical wave with an amplitude that decays steeply with increasing distance from the trapped bead (Fig. 2F) and a phase shift that gradually increases with distance (Fig. 2D).

To explain these results, we assume that the colloidal crystal can be treated as a 2D continuous elastic material. We write an equation of motion for the displacement field

$\overrightarrow{\mathbf{u}}$ in the solid (*18*) to which we add a dissipative term to account for the damping fluid and an oscillating point force,

, that represents the perturbation

$$\mathrm{\rho}\frac{{\mathrm{\partial}}^{2}\overrightarrow{\mathbf{u}}}{\mathrm{\partial}{t}^{2}}+\mathrm{\gamma}\frac{\mathrm{\partial}\overrightarrow{\mathbf{u}}}{\mathrm{\partial}t}=\overrightarrow{\mathbf{f}}(t)\mathrm{\delta}(\overrightarrow{r})+\frac{E}{2(1+\mathrm{\nu})}{\overrightarrow{\nabla}}^{2}\overrightarrow{\mathbf{u}}+\frac{E}{2(1-\mathrm{\nu})}\overrightarrow{\nabla}(\overrightarrow{\nabla}\xb7\overrightarrow{\mathbf{u}})$$(1)

Here, the first term describes the inertial forces with ρ ≈ 10^{−2} kg/m^{2} as the density of the 2D material. The second term represents the viscous damping due to the solvent with γ as the drag coefficient per unit area, which we determine experimentally from the short-time diffusion of particles in the crystal, giving γ = 4.8 · 10^{3} Ns/m^{3} (section S3). The last two terms correspond to the Navier-Cauchy equation that describes the elastic forces within the 2D solid with *E* as the 2D elastic modulus and ν as the 2D Poisson ratio.

Because our experiments are performed at low Reynolds number, the inertial term is negligible, resulting in overdamped mechanics. Solving the equation of motion for this case yields the displacement field in the form

$\overrightarrow{\mathbf{u}}=\mathbf{\alpha}\xb7\overrightarrow{\mathbf{f}}$, with **α** as the complex response function tensor, which has principal components α_{∥} and α_{⊥} that describe the components of the displacement field parallel and perpendicular to the applied force, respectively. In polar coordinates, with *r* as the distance from the point where the force is applied and θ = 0 corresponding to the direction of the force, this becomes (see section S4 for full details)

(2)and

$${\mathrm{\alpha}}_{\perp}(r,\mathrm{\theta})=\frac{1-{\mathrm{\nu}}^{2}}{4\mathrm{\pi}E}\text{sin}(2\mathrm{\theta})\left[{K}_{2}\right(\frac{r\sqrt{i}}{\mathrm{\zeta}})-{\mathrm{\lambda}}^{2}{K}_{2}(\frac{r\mathrm{\lambda}\sqrt{i}}{\mathrm{\zeta}}\left)\right]$$(3)where *K*_{0} and *K*_{2} denote the modified Bessel functions of the second kind,

is a characteristic attenuation length of the displacement amplitude, and

$\mathrm{\lambda}=\sqrt{2/(1-\mathrm{\nu})}$is a parameter that depends only on the Poisson ratio and diverges for ν → 1, which is the maximum Poisson ratio in 2D. The amplitude and the phase of the displacement fields are obtained as the magnitude and the argument, respectively, of the complex response functions.

To compare our experimental results with this prediction, we decompose the measured displacement amplitudes into their parallel and perpendicular components (Fig. 3, A and E). This results in distinct lobed patterns for these two components that can be observed in all our experiments. Our theoretical prediction produces identical patterns that are in excellent agreement, indicating that wave propagation in these colloidal crystals can indeed be described by treating the material as a continuous and isotropic 2D elastic solid (Fig. 3, B and F). We note that, in principle, the mechanical properties of the crystal should depend on the direction with respect to the crystal axes, but apparently, such anisotropy does not play a large role for these soft crystals.

The parallel component of the displacement response propagates preferentially along the excitation axis and shows a distinct asymmetry in the attenuation length along the two primary axes. The perpendicular displacement shows a four-lobed pattern, with maximum displacements at an angle of 45^{∘} with respect to the excitation direction. Also, in the phase maps (Fig. 3, C and D), we observe patterns that are in excellent agreement with the theoretical prediction (Fig. 3, G and H).

We can now use our theoretical analysis to interpret the experimental results in terms of the linear elasticity of the solid. For this, we first consider the parallel displacement component in the direction of the excitation, θ = 0. An asymptotic expansion of Eq. 2 for relatively large distances from the perturbation *r* > ζ (see section S4) leads to a phase lag in the far field

(4)and an amplitude

$${A}_{\parallel}(r,0)\sim {r}^{-1/2}{e}^{-\mathrm{r}/\mathrm{\zeta}\sqrt{2}}$$(5)

According to Eq. 4, the phase varies linearly with *r* along θ = 0, which is indeed what we find experimentally (Fig. 4A). This means that the wave propagates at a constant velocity in the direction of the excitation, with a phase velocity

. As shown in the inset of Fig. 4A, the phase velocity increases approximately as

${v}_{\parallel}\sim \sqrt{\mathrm{\omega}}$for low frequencies, which indicates that the elastic modulus and the Poisson ratio do not depend on the frequency in this regime.

The linear elasticity of an isotropic 2D solid is described by two independent mechanical parameters: the elastic modulus and the Poisson ratio. While the modulus only affects the absolute magnitude of the response function and the characteristic length scale ζ, the shape of the spatial pattern is uniquely determined by the Poisson ratio, as expressed by the parameter λ in Eqs. 2 and 3 (section S5). Hence, if the Poisson ratio is independent of the frequency, it should be possible to superimpose the measured displacement data obtained at different frequencies by normalizing the distance *r* by the characteristic length ζ. We determine ζ for each frequency from the asymptotic behavior at large *r* in two independent ways, using Eqs. 4 and 5. Normalizing all distances with ζ(ω) then indeed leads to a collapse of the data for both the amplitude and phase of the parallel displacement contribution, both in the direction of the excitation (θ = 0) and perpendicular to it (θ = π/2) (Fig. 4, B and C). We fit these curves to Eq, 2 using

as the only fit parameter, giving a value for the Poisson ratio ν = 0.70 ± 0.07, comparable to values expected for crystalline solids in two dimensions (*19*, *20*). Using this value of ν, we then obtain the elastic moduli of the colloidal crystal from the decay lengths, giving moduli on the order of *E* ≈ 10^{−6} N/m, both for the phase and amplitude data (Fig. 4D). A similar analysis for the perpendicular displacement components, fitted to Eq. 3 (section S6), gives similar moduli (Fig. 4D). As a final consistency check, we use the measured Young’s modulus and Poisson ratio to predict the absolute values of the response functions, which now provide a reasonable quantitative match with the experimental results (Fig. 3, A to H). This highlights that our theoretical description captures the main phenomena in a quantitative fashion.

We may compare the values for *E* with a simple estimate obtained by approximating the colloidal crystal as a hexagonal lattice of harmonic springs, for which

with *k* as the spring constant of a particle pair (*21*). We estimate *k* by analyzing the thermal bond length fluctuations (section S7) and find *E* ≈ 1 · 10^{−6} N/m in very good agreement with the experimental values. We also validate our results by estimating the modulus in a different way, by analyzing the long-wavelength thermal fluctuations in the crystal (*22*). As shown in section S8, this gives a modulus that is in good agreement with our estimate based on the wave pattern. Our method has the advantage, however, that it can directly probe the frequency dependence of the elastic moduli. Moreover, we obtain a spatial map of the displacement pattern in real space, which allows us to observe deviations from the predicted pattern and to relate these to the local structure (see below).

We note that our method is limited at higher frequencies by the decrease in the attenuation length with increasing frequency. Once the characteristic length ζ becomes of the order of the particle size, discretization effects hinder the accurate determination of the wave propagation. This gives a limiting frequency ω_{max} ≈ *E*/γ*a*^{2} ≈ 20 rad/s. For frequencies approaching ω_{max}, our estimation of the modulus becomes inaccurate.

Last, to emphasize that our approach is not exclusive to ordered, crystalline solids, we repeat the experiments described above for a disordered colloidal glass that lacks long-ranged order (*23*, *24*). We prepare monolayers of a bidisperse mixture of silica spheres (*a* = 3.04 and 1.85 μm at a mixing ratio of 1:2.35 by number and a total packing fraction of 0.80) (Fig. 1B). We confirm the absence of structural order from the liquid-like shape of the pair-correlation function and structure factor, while the particle dynamics are strongly arrested and caged as evidenced from their mean square displacement in the absence of external mechanical excitation (section S1.2).

Despite the very different microstructure, we find that the displacement amplitude and phase patterns are very similar to those observed for the crystals (Fig. 5, A to D). The same analysis (Fig. 5, E and F) as for the crystals (using a drag coefficient γ that is a weighted average over the small and large particles, giving γ = 7.5 · 10^{2} Ns/m^{3}; see section S3) gives a Poisson ratio ν = 0.60 ± 0.15 and an elastic modulus on the order of 1 · 10^{−7} N/m for the glasses (log *E* = − 7.0 ± 0.4). This is about an order of magnitude lower than for the crystals, which could be due to the slightly lower packing fraction or the softer interaction potential between the small particles (section S7).

From Fig. 5 (A to D), it is also clear that the wave patterns are noisier for the glasses than for the crystals, even though the decay lengths ζ are comparable. The reason for this is not clear, but it might be related to their disordered structure, which is known to lead to inhomogeneous mechanical properties and nonaffine deformation modes (*25*). Further evidence for a relation between wave propagation and microscopic structure can be seen by looking at the displacement patterns in the vicinity of defects in the colloidal crystal. For example, the pattern in Fig. 2F shows an asymmetry between left and right, which can be attributed to the presence of a dislocation just to the right of the probe particle (section S9). We note that this is not yet captured by the theory, which assumes a homogeneous elasticity.

Our results highlight that the observation of highly damped mechanical waves at low *Re* open up the way for mechanical characterization of both ordered and disordered ultraweak solids. Traditionally, microrheology is the method of choice to characterize the viscoelasticity of very weak elastic materials whose moduli are below the detection limit of conventional macroscopic rheometers. In microrheology, the viscoelastic features of the material are extracted from the Brownian fluctuations of tracer particles embedded in the material (passive microrheology) (*26*) or from the response of an actively forced probe particle (active microrheology) (*27*). However, these methods often fail to give the bulk rheological properties, for example, in very heterogeneous samples where the probe particles sample mostly the pores in the material or in cases where the probe particles locally modify the structure of the material. This can be solved by analyzing cross-correlations between particles that are separated by a distance much larger than the characteristic length scale of the material, in microrheology is a field of science so-called two-particle microrheology (*28*). Alternatively, the elastic properties can be obtained by analyzing the long-wavelength thermal vibrational modes of the material using particle tracking techniques, although extracting the frequency-dependent rheological properties is more challenging in this case (*22*). Our method can be seen as a form of active multiparticle microrheology, in which the complete spatial displacement pattern induced by a local perturbation is mapped out. By filtering out the Brownian noise, we obtain experimental access to the frequency-dependent mechanical properties of extremely soft and fragile systems. In principle, the method can also be extended toward the nonlinear deformation regime by increasing the amplitude of the probe oscillation, although the interpretation in terms of nonlinear elastic moduli will obviously be more challenging.

**Acknowledgments: ****Funding:** This work is part of the Industrial Partnership Program Hybrid Soft Materials that is carried out under an agreement between Unilever Research and Development B.V. and the Netherlands Organisation for Scientific Research (NWO). J.v.d.G. acknowledges the European Research Council for financial support (ERC Consolidator grant Softbreak). The work of J.S. is part of the VIDI research program (no. 723.016.001) of NWO. **Author contributions:** J.M.v.D. and J.S. conceived the experiments. J.M.v.D., R.W., and R.F. constructed and programmed the optical tweezer instrument. J.M.v.D. and R.H. developed the real-time distributed particle tracking protocol. J.M.v.D. performed the experiments and data analysis. J.v.d.G. developed the wave propagation theory. J.M.v.D., J.v.d.G., A.Z., and J.S. discussed the data interpretation and cowrote the manuscript. **Competing interests:** The authors declare that they have no competing interests. **Data and materials availability:** All data needed to evaluate the conclusions in this paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.