Fast spectral solver for viscoelastic structures under oscillatory flow in free space or wall-bounded domains: applications to quartz crystal microbalance and force spectroscopy

TL;DR

This work presents a fast spectral solver for the linear viscoelastic response of immersed structures in oscillatory flow, formulated in the frequency domain to connect viscoelastic networks with oscillatory Stokes flow in open, doubly periodic domains. It offers two coupling strategies—the mobility-based route and an Anderson-accelerated iterative scheme—implemented via a highly efficient, space- and time-spectral fluid solver. The approach is validated through quartz crystal microbalance (QCM) scenarios with rigid and viscoelastic linkers, and force spectroscopy by AFM, achieving quantitative agreement with analytical, numerical, and experimental benchmarks across scales from nanometers to microns. The framework enables quantitative inference of intrinsic molecular viscoelastic properties from experimental signals by bridging mesoscopic models and first-principles descriptions, with significant potential for integrated theory–experiment analyses in soft matter and biophysics.

Abstract

We present a fast spectral solver for the linear response of viscoelastic structures under oscillatory flow either in free space or close to a flat moving wall. The scheme works in the frequency domain (using phasors) and couples the oscillatory Stokes equation with rigid or flexible structures, modeled by viscoelastic networks of immersed boundary kernels. The fluid-structure coupling can be solved by two routes. One route calculates the hydrodynamic mobility matrix required to solve the equation for the structure deformation rate in matrix form. The second route iteratively solves the coupled fluid-structure equations: fluid-induced forces on the structures create a tension field which is then transferred to the fluid, until convergence. The resulting fixed-point problem is solved iteratively using the Anderson acceleration method. The mobility route is optimal when dealing with one or few structures, while the iterative scheme is preferred for denser dispersions. In any case, the flow resulting from the body forces is solved by a recently developed scheme [J. Fluid. Mech. 1010 A57, 2025] which is spectral in space and time and deals with doubly periodic open domains (either free-space or wall-bounded) where meshing is restricted to the region of interest around the structures. We test the present scheme in two applied contexts: quartz-crystal-microbalance (QCM) of spheres, suspended, adsorbed or tethered to viscoelastic linkers; and force spectroscopy (via atomic force microscopy) reproducing the power spectra of vibrating microparticles near a solid boundary. In all cases, comparisons with analytical, numerical and experimental results show excellent agreement. We conclude by discussing new routes the scheme opens in force spectroscopy and QCM analyses of soft objects.

Figures (13)

• Figure 1: Two types of experiments reproduced by the present solver. (Left) Quarzt crystal microbalance experiments and (right) vibrational spectra obtained from atomic force microscopy (AFM). The left panel schematically represent three different levels of description for the same system: a FtsZ protein having a flexible peptide chain along its intrinsically disordered region (IDR). (Top-left) The microscale, all-atom model (water molecules not shown and colors indicate partial charges of the residues) can provide information to the mesoscale model (center) by a Bottom-Up (BU) formalism. The mesoscale level is based on a coarse-grained (CG) model with a globular part (radius $R_g$) formed by a multiblob structure (12 blobs in the vertices and one in the center of an isoshedron). The later is connected to the highly flexible IDR, modelled by a viscoelastic spring . The CG model is coupled to the oscillatory solvent flow via immersed-boundary continuum fluid dynamics 2025_spectral_solver. The CG parameters are the intrisinc friction of the chain connected to the substrate $\xi$, the chain stiffness $k$ and average length $L$, the chain orientation angles, as well as the elastic network of springs between the beads forming the rigid globule. At the macroscale (top-right corresponding to the QCM experiment) a Top-Down (TD) approach provides connection with the meso-level, via quantitative reproduction of experimental signals (frequency and dissipation shifts, directly related to the complex-phasor of the averaged excess surface tangential stress $\bar{\sigma}_{xz}$). In the QCM device, oscillations of the quartz crystal substrate are excited by an oscillatory voltage, via the inverse-piezo-electric effect, leading to the substrate horizontal oscillation at angular frequency $\omega_n \in 10 \pi \, n \, 10^{6}\, \mathrm{rad/s}$ and nanometer amplitude ($n=\left\{1,3,5,..\right\}$ is called the overtone)
• Figure 2: Real (red) and imaginary (blue) components of the impedance at the suspended particle wall as a function of height, calculated using the method described in this work (dots) and compared to the results of Leshansky (2023)Leshansky2023. For the numerical calculations, the sphere, of geometric radius $R_g$ (see appendix \ref{['app:Lagrange']}) is represented using 92 3-point Peskin kernels. The box side was approximately $L_\parallel\approx8\delta$, and the meshed domain $z\in [-H,H]$ in the unbounded $z$ direction was adapted to the maximum height in each case as $H=\max(z_0)+2R_g$. The impedance is scaled with $\eta\, \tilde{n} R_g$ where $\tilde{n}$ is the particle surface number density.
• Figure 3: Real (red) and imaginary (blue) components of the impedance at the wall for adsorbed rigid particles of geometrical radius $R_g$, calculated using the method described in this work (dots), compared to the theoretical predictions from Ref. Leshansky2024. The number of blobs per sphere is 1415, and the location of the IB blobs follows the protocol explained in the appendix \ref{['app:Lagrange']} for filled spheres with $T=6$.
• Figure 4: Real (left panel) and imaginary (right panel) components of the wall impedance for a particle connected to the wall via a viscoelastic linker, shown as a function of the damping coefficient for different values of the elastic constant. Results from the numerical code (dots) are compared with the predictions of an analytical theory (lines). Calculations were performed for neutrally buoyant blobs with a hydrodynamic radius $a = 0.03\delta$, a linker length $l_{\text{link}} = 0.15\delta$, with an longitude and azimuth angle of $45^\circ$. The box size is $L = 50a$, ensuring that the effects of periodic boundary conditions are negligible.
• Figure 5: Contour plot for the frequency shift $\Delta F/n$ (Hz) and the disipation $\Delta D\times10^6$ scaled with the coverage $\Theta$, for a viscoelastic linker connected to a bead of $b=2.5\text{nm}$ radius. The linker model corresponds to a Gaussian linker with average lenght $\langle \ell \rangle =10\text{nm}$ and standard deviation $\sigma_\ell= 3\,\text{nm}$. The x-axis corresponds to the elasticity parameter of the linker $\kappa=k/(\omega \xi_S)$ and the y-axis to the internal friction $\gamma=\xi/\xi_S$, the bead mass is $m^e=0.2m_f$ (compatible with a protein). The impedance is averaged according to eq. \ref{['eq:zave']} for the $n=7$ overtone of a QCM resonator with fundamental frequency $f_0=5\text{MHz}$ in water ($\nu=10^{-6}m^2/s$). The red rectangles indicate the region of parameters for a flexible peptide chain with 60 aminoacids similar to the disordered domain of the FtsZ protein illustrated in Fig. \ref{['fig:levels']}(a): contour length $L_c=20 \text{nm}$, persistence length $0.35 \text{nm}$ and stiffness $k\approx 0.4 k_B/\text{nm}^2$, which corresponds to $\kappa \approx 1$ for the $n=7$ overtone.