Efficient numerical methods for the Maxey-Riley-Gatignol equations with Basset history term

TL;DR

The Maxey-Riley-Gatignol equations with the Basset history term pose significant memory and nonlocal challenges. The paper adopts a diffusion-type reformulation on a semi-infinite domain and develops second- and fourth-order finite-difference schemes to solve it, paired with DIRK and IMEX Runge-Kutta time stepping to avoid nonlinear solves. A comprehensive comparison against Daitche’s direct integrator and Prasath’s polynomial-expansion method shows that order of convergence can degrade for nonzero initial relative velocities or nonneutral buoyancy, but IMEX strategies offer real-time capability in several regimes, and the FD4 scheme excels for zero relative velocity. The study demonstrates that solving the full MRGE in real time is feasible for certain flow fields (e.g., Bickley jet, Faraday) and particle properties, highlighting practical pathways for particle-tracking applications in complex flows.

Abstract

The Maxey-Riley-Gatignol equations (MRGE) describe the motion of a finite-sized, spherical particle in a fluid. Because of wake effects, the force acting on a particle depends on its past trajectory. This is modelled by an integral term in the MRGE, also called Basset force, that makes its numerical solution challenging and memory intensive. A recent approach proposed by Prasath et al. exploits connections between the integral term and fractional derivatives to reformulate the MRGE as a time-dependent partial differential equation on a semi-infinite pseudo-space. They also propose a numerical algorithm based on polynomial expansions. This paper develops a numerical approach based on finite difference instead, by adopting techniques by Koleva et al. and Fazio et al. to cope with the issues of having an unbounded spatial domain. We compare convergence order and computational efficiency for particles of varying size and density of the polynomial expansion by Prasath et al., our finite difference schemes and a direct integrator for the MRGE based on multi-step methods proposed by Daitche. While all methods achieve their theoretical convergence order for neutrally buoyant particles with zero initial relative velocity, they suffer from various degrees of order reduction if the initial relative velocity is non-zero or the particle has a different density than the fluid.

Figures (11)

• Figure 1: Trajectories of two particles moving in the vortex flow field with initial position $\bm{y}(t_0) = (1,0)^T$ and zero initial relative velocity with $S=0.3$ and $R=7/9$ (left) or $R=4/3$ (right) for $t\in[0,10]$.
• Figure 2: Trajectories of two particles moving in the oscillatory flow field with initial position $\bm{y}(t_0) = {(0,0)}^T$ and zero initial relative velocity with $S=0.3$ and $R=7/9$ (left) or $R=4/3$ (right) for $t\in[0,3]$. The arrows show the velocity field at $t=3$.
• Figure 3: Horizontal component of the position versus time for two particles decelerating in a quiescent flow with initial position $\bm{y}(t_0) = {(0,0)}^T$ and initial relative velocity $\bm{q}(0, t_0) = {(0.1, 0)}^T$ with $S=0.3$ and $R=7/9$ (left) or $R=4/3$ (right) for $t\in[0,1]$.
• Figure 4: Trajectory of two particles moving in the Bickley jet with initial position $\bm{y}(t_0) = {(0, 0)}^T$ and zero initial relative velocity with $S=0.3$ and $R=7/9$ (left) and $R=4/3$ (right) for $t\in[0,1]$. The arrows show the velocity field at $t=1$.
• Figure 5: Trajectory of two particles moving in the Faraday flow with initial position $\bm{y}(t_0) = {(0.02, 0.01)}^T$ and zero initial relative velocity, i.e. $\bm{q}(0, t_0) = {(0, 0)}^T$, with $S=0.3$ and either $R=7/9$ (left) or $R=4/3$ (right) for $t\in[0,5]$. The arrows show the velocity field at $t=5$.