Generating, from scratch, the near-field asymptotic forms of scalar resistance functions for two unequal rigid spheres in low Reynolds number flow

The motion of rigid spherical particles suspended in a low Reynolds number fluid can be related to the forces, torques and stresslets acting upon them by 22 scalar resistance functions, commonly notated $X^A_{11}$, $X^A_{12}$, $Y^A_{11}$, etc. Near-field asymptotic forms of these resistance functions were derived in Jeffrey and Onishi (J. Fluid Mech., 1984) and Jeffrey (Phys. Fluids A, 1992); these forms are now used in several numerical methods for suspension mechanics. However, the first of these important papers contains a number of small errors which make it difficult for the reader to correctly evaluate the functions for parameters not explicitly tabulated. This short article comprehensively corrects these errors, and adds formulae that were originally omitted from both papers, so that the reader can verify and implement the equations independently. The corrected expressions, rationalised and using contemporary nondimensionalisation, are shown to match mid-field values of these scalars which are calculated through an alternative method. A Python script to generate and evaluate these functions is provided.


Introduction
The linearity of Stokes flow allows the forces, torques and stresslets acting on particles suspended in a low Reynolds number fluid to be related to the particles' velocities, angular velocities and fluid background strain rate through a resistance or mobility matrix.For rigid spherical particles, the elements of the resistance matrix can be formed from scalar resistance functions, commonly notated as The near-field asymptotic forms of the scalar resistance functions for two unequal rigid spheres were initially published in Jeffrey & Onishi (1984) and Jeffrey (1992).The terms in the asymptotic forms themselves rely on recurrence relations, and, in the first of these two papers, there are a number of small errors.Furthermore, the reader may find it difficult to generate the required expressions (and therefore the value of these functions) in both papers independently, given the omission of intermediate formulae.
Some of these errors appear to have been noticed by authors, including indirectly by Kim & Karrila (2005) in corrected tabulations, and a partial list of errata has been published online by Ichiki (2008).However, a comprehensive, fully corrected description of how to generate expressions for these functions, to enable verifiable calculation for arbitrary size ratios and separation distances, is yet to be published.This short article is a compilation of the relevant equations, with those originally omitted now added, and with any errors fixed.Equations from Jeffrey & Onishi (1984) are labelled as (JO 1.1); similarly, equations from Jeffrey (1992) are labelled as (J 1), and those from the helpful Ichiki et al. (2013) are labelled as (I 1).Pages from Kim & Karrila (2005) are labelled as (K&K, p. 1).
At the end of the article, the questions of error identification and validation of the corrected forms are addressed.
Throughout this article, we use the same notation as the original papers.For two spheres of radius a 1 , a 2 with centres a distance s apart, we define the non-dimensional gap ξ and size ratio λ as The near-field forms are valid for ξ ≪ 1 and ξ ≪ λ.
The nondimensionalisation in (JO 1.7) and (J 3) scales the resistance functions by factors of (a 1 + a 2 )/2.As noted in Jeffrey (1992), although this choice provides symmetry elsewhere, it leads to unwieldy factors of (1 + λ) −1 in some of the expressions we are interested in.Here, we instead follow the convention of Kim & Karrila (2005, §11.3), and scale all functions by factors of a 1 .Specifically, to recover the dimensional forms, it is necessary to multiply the scalars in eq. ( 1), derived as follows, by an appropriate factor: A Python script implementing the corrected formulae is available on GitHub (Townsend, 2023).

X A terms
Here the X A formulae are given in full, with changes from the source material when noted.
The same directions for alteration, when required, will be given for the other terms in later sections.
We first set up up the recurrence relations from (JO 3.6-3.9), V npq = P npq − 2n (n + 1)(2n + 3) Then, we define the formulae from (JO 3.15 and 3.19), We then define two intermediate functions, the first of which is given in (JO 3.22) as The function m 1 is somewhat awkward as, despite being a toggle for whether m = 2, it is used in Jeffrey & Onishi (1984) for summation over both even and odd values of m.The later paper, Jeffrey (1992), removes it by using instead a form of the summation which here would render the equivalent expression It is computationally convenient to have the same coefficients of f m , g 1 , g 2 and g 3 (where they exist) in the summations for all functions of this type, so we choose to use this form when describing corrections going forward.Consistent coefficients also allow us to notice that all f m functions are defined with a factor of 2 m which cancels the 2 −m in the f m coefficient in the sum: these powers of 2 can therefore be left out of numerical implementations.In practice, we also find that moving the factor (1 + λ) −m inside the sum in the definition of f m can prevent the terms in the sum (specifically from λ q ) from exceeding the maximum float size when λ and m (and hence q) are large.However, we do neither of these things in this paper to remain consistent with the original texts.
The second intermediate function is from (JO 3.23), and indeed needs correcting to be noting the change from m 1 (which is identically m − 2 given that m is odd) to m + 2.
Then the resistance scalars are given by from (JO 3.17-3.18)up to O(ξ log(ξ −1 )): note the different scaling on X A 12 as we use the aforementioned convention from (K&K, p. 279).

Y A terms
The recurrence relations are , but with V npq corrected to noticing the sign change on the 1 in the last subscript.
Then the A Y terms are given by (JO 4.17-4.18):if desired, these can be more conveniently written to remove m 1 and use the Jeffrey (1992) form by writing but this does not constitute a correction.
Either way, we are left with the resistance scalar formulae, from , with a different scaling on Y A 12 .

Y B terms
The recurrence relations are the same as those for the Y A terms.The required intermediate formulae for the f and g functions are and (JO between 5.6 and 5.7), respectively.
The B Y terms are corrected from (JO 5.7-5.8) to become Then the resistance scalars are given by from (JO 5.5-5.6), with a different scaling on Y B 12 .

X C terms
Expressions for the resistance scalars can be expressed directly as from (JO 6.9-6.10),where X (28)

Y C terms
The recurrence relations are the same as those for the Y A terms except that the initial conditions are replaced by ).
The intermediate formula for the f function is with the g formula given by (JO between 7.10 and 7.11), with the correction to g 5 of Then the C Y terms from , correcting the latter [note the power (1 + λ)] and writing both to avoid m 1 , are C 12 uses the different scaling, and where ζ(z, a) is the Hurwitz zeta function,