Karthik Sasihithlu
email: karthik.sasihithlu@institutoptique.fr

Arvind Narayanaswamy

Near-field radiative transfer between two unequal sized spheres with large size disparities

We compute near-field radiative transfer between two spheres of unequal radii R 1 and R 2 such that R 2 40R 1 . For R 2 = 40R 1 , the smallest gap to which we have been able to compute radiative transfer is d = 0.016R 1 . To accomplish these computations, we have had to modify existing methods for computing near-field radiative transfer between two spheres in the following ways: (1) exact calculations of coefficients of vector translation theorem are replaced by approximations valid for the limit d R 1 , and (2) recursion relations for a normalized form of translation coefficients are derived which enable us to replace computations of spherical Bessel and Hankel functions by computations of ratios of spherical Bessel or spherical Hankel functions. The results are then compared with the predictions of the modified proximity approximation.

Introduction

Theoretical and experimental investigation of near-field radiative enhancement between bodies has been an area of much study recently due to its possible application in a range of areas including thermal photovoltaics [START_REF] Narayanaswamy | Surface modes for near field thermophotovoltaics[END_REF][START_REF] Basu | Review of near-field thermal radiation and its application to energy conversion[END_REF][START_REF] Ilic | Overcoming the black body limit in plasmonic and graphene near-field thermophotovoltaic systems[END_REF][START_REF] Laroche | Near-field thermophotovoltaic energy conversion[END_REF], thermal rectification [START_REF] Otey | Thermal rectification through vacuum[END_REF], nanopatterning [START_REF] Lee | Confinement of infrared radiation to nanometer scales through metallic slit arrays[END_REF], thermally assisted magnetic recording [START_REF] Challener | Heat-assisted magnetic recording by a near-field transducer with efficient optical energy transfer[END_REF], thermal imaging [START_REF] Wilde | Thermal radiation scanning tunnelling microscopy[END_REF] and noncontact radiative cooling [START_REF] Guha | Near-field radiative cooling of nanostructuresm[END_REF]. The theoretical analysis has relied on Rytov's theory of fluctuational electrodynamics [START_REF] Rytov | Theory of electric fluctuations and thermal radiation[END_REF] and has been used to compute near-field radiative transfer between two flat surfaces [START_REF] Polder | Theory of radiative heat transfer between closely spaced bodies[END_REF][START_REF] Loomis | Theory of heat transfer by evanescent electromagnetic waves[END_REF][START_REF] Cravalho | Effect of small spacings on radiative transfer between two dielectrics[END_REF][START_REF] Mulet | Enhanced radiative transfer at nanometric distances[END_REF][START_REF] Pendry | Radiative exchange of heat between nanostructures[END_REF][START_REF] Park | Proximity-enhanced thermal radiation[END_REF], between two dipoles [START_REF] Volokitin | Radiative heat transfer between nanostructures[END_REF][START_REF] Domingues | Heat transfer between two nanoparticles through near-field interaction[END_REF], and between a dipole and a flat surface [START_REF] Pendry | Radiative exchange of heat between nanostructures[END_REF][START_REF] Mulet | Nanoscale radiative heat transfer between a small particle and a plane surface[END_REF]. Only recently has this been applied to other more experimentally feasible geometries such as between two mesoscopic spheres [START_REF] Narayanaswamy | Thermal near-field radiative transfer between two spheres[END_REF] and between a mesoscopic sphere and a planar surface [START_REF] Otey | Numerically exact calculation of electromagnetic heat transfer between a dielectric sphere and plate[END_REF][START_REF] Krüger | Nonequilibrium electromagnetic fluctuations: Heat transfer and interactions[END_REF]. The adjective mesoscopic is being used here to describe objects whose sizes are of the order of the characteristic thermal wavelength of radiation λ T (λ T ≈ hc/k B T , where 2π h is the Planck's constant, c is the speed of light in vacuum, k B is the Boltzmann constant, and T is the absolute temperature). In contrast, the size of macroscopic objects is much greater than λ T . The formalism for computing radiative transfer between arbitrary bodies has also been detailed in several recent publications [START_REF] Krüger | Nonequilibrium electromagnetic fluctuations: Heat transfer and interactions[END_REF][START_REF] Rodriguez | Fluctuating-surface-current formulation of radiative heat transfer for arbitrary geometries[END_REF][START_REF] Messina | Scattering-matrix approach to Casimir-Lifshitz force and heat transfer out of thermal equilibrium between arbitrary bodies[END_REF][START_REF] Narayanaswamy | A Green's function formalism of energy and momentum transfer in fluctuational electrodynamics[END_REF][START_REF] Ben-Abdallah | Many-body radiative heat transfer theory[END_REF]. Recent experimental validation of this near-field enhancement for sub-micron gaps has been limited to that between a mesoscopic sphere and a plate [START_REF] Shen | Surface phonon polaritons mediated energy transfer between nanoscale gaps[END_REF][START_REF] Rousseau | Radiative heat transfer at the nanoscale[END_REF][START_REF] Van Zwol | Tuning near field radiative heat flux through surface excitations with a metal insulator transition[END_REF][START_REF] Shen | Nanoscale thermal radiation between two gold surfaces[END_REF][START_REF] Narayanaswamy | Near-field radiative heat transfer between a sphere and a substrate[END_REF] and between a scanning probe tip and a planar surface [START_REF] Hirsch | Heat transfer in ultrahigh vacuum scanning thermal microscopy[END_REF][START_REF] Kittel | Near-field heat transfer in a scanning thermal microscope[END_REF][START_REF] Worbes | Enhanced near-field heat flow of a monolayer dielectric island[END_REF].

The exact computation of the near-field thermal radiative energy transfer between two adjacent bodies involves expressing the electric and magnetic dyadic Green's function (DGF) [START_REF] Yaghjian | Electric dyadic Green's functions in the source region[END_REF][START_REF] Collin | Field Theory of Guided Waves[END_REF] of the vector Helmholtz equation as a sum of properly weighted vector eigenfunctions (e.g. vector spherical waves in spherical coordinates or vector planar waves in Cartesian coordinates) of the vector Helmholtz equation. The weights of the vector eigenfunctions are determined by imposing the appropriate boundary conditions on the surfaces of the bodies. The boundary conditions account for multiple reflections and scattering of the outgoing thermal radiation from the bodies. The relation between the thermally generated electric currents, which are the source for thermal radiation, and the temperature of the bodies is given by the fluctuation-dissipation theorem [START_REF] Landau | Electrodynamics of Continuous Media[END_REF][START_REF] Kubo | The fluctuation-dissipation theorem[END_REF][START_REF] Marconi | Fluctuation-dissipation: response theory in statistical physics[END_REF]. This procedure has been adopted to compute the thermal near-field radiative transfer between two spheres [START_REF] Narayanaswamy | Thermal near-field radiative transfer between two spheres[END_REF] and between a sphere and a plane [START_REF] Otey | Numerically exact calculation of electromagnetic heat transfer between a dielectric sphere and plate[END_REF].

Crucial to the analysis of radiative transfer between two spheres is the computation of coefficients of the translation addition theorem to facilitate the application of boundary conditions on the surface of the adjacent objects. In particular, the vector addition theorem allows us to transform vector spherical waves in a coordinate system with origin at the center of one of the spheres into vector spherical waves in a coordinate system with origin at the center of the other sphere. While there exist exact expressions to compute the coefficients of the vector addition theorem [START_REF] Bruning | Multiple scattering by spheres[END_REF][START_REF] Chew | Derivation of the vector addition theorem[END_REF], recursion relations have to be used to ensure their efficient computation [START_REF] Chew | Efficient ways to compute the vector addition theorem[END_REF].

In addition to the exact method [START_REF] Narayanaswamy | Thermal near-field radiative transfer between two spheres[END_REF] for computing radiative transfer between spherical objects, which can be computationally expensive at small gaps, there exists an approximate method to predict the radiative transfer between nonplanar surfaces by treating them as a series of parallel surfaces with varying gaps. This method, called the modified proximity approximation method (MPA) [START_REF] Sasihithlu | Proximity effects in radiative heat transfer[END_REF] differs from the proximity force approximation (PFA) [START_REF] Błocki | Proximity forces[END_REF][START_REF] Lamoreaux | Demonstration of the Casimir force in the 0.6 to 6µm range[END_REF] used in Casimir physics in only the way the far-field contribution to the radiative transfer (which is negligible for forces) is treated. While there has been no rigorous proof in literature in deriving such a formulation, the validity of the PFA in different limits has been investigated recently for Casimir force [START_REF] Gies | Casimir effect for curved geometries: Proximity-force-approximation validity lim-its[END_REF] and the near-field radiative transfer [START_REF] Golyk | Small distance expansion for radiative heat transfer between curved objects[END_REF] between a sphere and a plane. The MPA was first developed to ensure continuity in the radiative transfer between two equal sized spheres in the far-field and the near-field limits [START_REF] Sasihithlu | Proximity effects in radiative heat transfer[END_REF]. Hence it would be of interest to also study the relevance of the MPA when the radius of one of the spheres R 2 is increased such that R 2 R 1 where R 1 is the radius of the smaller sphere, since when R 2 /R 1 → ∞ it reduces to the sphere-plane configuration that has been adopted for experimental verification.

Though the theoretical framework for the analysis of near-field radiative transfer between two spheres has been established by Narayanaswamy and Chen [START_REF] Narayanaswamy | Thermal near-field radiative transfer between two spheres[END_REF], it cannot be extended as is for unequal sized spheres with large radius ratio due to computational constraints. It is our intention in this work to develop simplifications to the computational model which would enable us to extend it to the case of disproportionately sized spheres. The paper has been organized as follows. In Sec. 2 we give a brief introduction to the translation addition theorem since it is central to our work and list the recurrence relations that have been developed in literature to compute the translation coefficients. We also propose a method to normalize the translation coefficients and show that the recursion relations for the computation of these normalized translation coefficients can be expressed in terms of ratios of spherical Bessel and Hankel functions. In Sec. 3, we propose a simplification for the translation addition theorem, which we have referred to as the one-term approximation, that aids the computation of near-field conductance between spherical bodies greatly. In Sec. 4, a description of the dependence of the normalized translation coefficients on the relative size of the two spheres is given which enables further simplification of the the computational procedure. In Sec. 5 we employ the simplifications listed above to compute near-field radiative transfer between two unequal sized spheres and compare this with the predictions from the MPA. We wish to stress that the approximations/methods developed here are have relevance not only to the computation of near-field radiative transfer but also to other problems in electromagnetic scattering, especially ones where electromagnetic surface waves and/or evanescent waves play an important role.

Translation addition theorems

Using the translation addition theorem, it is possible to transform the multipole expansion of an electromagnetic field from a coordinate system at the center of one scatterer to that of another in a multi-particle system. This becomes necessary while imposing the appropriate boundary conditions on the surface of the particles. Efficient computation of translation coefficients is essential, especially when the number of terms required for convergence in the multipole expansion is large, as is the case in near-field electromagnetic scattering problems. The simplest case of a multi-particle system is that of two spheres whose centers are translated along the zaxis as shown in Fig. 1. Emission and scattering of electromagnetic waves from such a system involves solving the vector Helmholtz equation given by: ∇ ∇ ∇ × ∇ ∇ ∇ × P P P( ( (r r r) ) )k 2 P P P( ( (r r r) ) ) = 0, (

where P P P( ( (r r r) ) ) is the electric or magnetic field at position vector r r r and k is the frequency dependent wave propagation constant. Following the notation in Fig. 1, k = k f = ω/c in the free-space region; k = k a = n a (ω) ω/c inside sphere 'a'; and k = k b = n b (ω) ω/c inside sphere 'b' with n a (ω) and n b (ω) being the complex refractive indices of sphere 'a' and sphere 'b' respectively. The independent divergence-free solutions of the vector Helmholtz equation in spherical coordinates are given by [START_REF] Narayanaswamy | Thermal near-field radiative transfer between two spheres[END_REF]:

M M M ( j) nm (kr r r) = z ( j) n (kr)V V V (2) nm (θ , φ ), (2)
N N N ( j) nm (kr r r) = ζ ( j) n (kr)V V V (3) nm (θ , φ ) + z ( j) n (kr) kr n(n + 1)V V V (1) nm (θ , φ ), (3)
where M M M ( j) nm (kr r r) and N N N

( j)
nm (kr r r) are vector spherical waves of order (n, m). n can take integer values from 0 to ∞. For each n, |m| ≤ n. The superscript p refers to the radial behavior of the waves. For j = 1, the M M M and N N N waves are regular waves and remain finite at the origin and z

j) n (x) = 1 x d dx xz ( j) n (x) . V V V (1) nm (θ , φ ), V V V (2) nm (θ , φ ), and V V V (3)
nm (θ , φ ) are vector spherical harmonics of order (n, m) and are given by:

V V V (1) nm (θ , φ ) = r r rY nm , (4a)
V V V (2) nm (θ , φ ) = 1 n(n + 1) -φ φ φ ∂Y nm ∂ θ + θ θ θ im sinθ Y nm , (4b)
and

V V V (3) nm (θ , φ ) = 1 n(n + 1) θ θ θ ∂Y nm ∂ θ + φ φ φ im sinθ Y nm . (4c)
Following Ref. [START_REF] Chew | Efficient ways to compute the vector addition theorem[END_REF] these vector spherical waves can also be expressed in terms of the solutions to the scalar Helmholtz equation. The solutions of the scalar Helmholtz equation ∇ 2 ψ + k 2 ψ = 0 are given by:

ψ ( j) nm (k, r r r) = z ( j) n (kr)Y nm (θ , φ ), (5)
where, ψ nm are scalar spherical waves of the same order (n, m) and Y nm (θ , φ ) represent the spherical harmonics. Using these the vector spherical waves can be expressed as :

M M M ( j) nm (k, r r r) = 1 n(n + 1) ∇ × r r rψ ( j) nm (k, r r r), (6)
and

N N N ( j) nm (k, r r r) = 1 n(n + 1) 1 k ∇ × ∇ × r r rψ ( j) nm (k, r r r). (7)
The representations for M M M and N N N waves in Eqs. ( 2), ( 3), [START_REF] Lee | Confinement of infrared radiation to nanometer scales through metallic slit arrays[END_REF], and (7) are equivalent. They differ from the corresponding representations in Ref. [START_REF] Chew | Efficient ways to compute the vector addition theorem[END_REF] in only a factor 1/ n(n + 1) which has

= r -R 1 -R 2 been introduced to ensure that V V V (1) nm (θ , φ ), V V V (2)
nm (θ , φ ), and V V V

nm (θ , φ ) are orthonormal to each other.

The scalar addition theorem for the scalar spherical waves ψ nm is given by:

ψ nm (k, r r r) = ∞ ∑ ν=1 ν ∑ m =-ν ψ νm (k, r r r )β νm ,nm , (8)
and the vector addition theorem for the vector spherical waves M M M ( j) nm and N N N

( j)
nm waves is given by:

M M M ( j) nm (kr r r) = ∞ ∑ ν=1 ν ∑ m =-ν M M M ( j) νm (kr r r )A νm ,nm + N N N ( j) νm (kr r r )B νm ,nm , (9)
N N N ( j) nm (kr r r) = ν=∞ ∑ ν=1 ν ∑ m =-ν N N N ( j) νm (kr r r )A νm ,nm + M M M ( j) νm (kr r r )B νm ,nm (10)
where r r r and r r r are coordinates with respect to the coordinate systems centered at O and O respectively as shown in Fig. 1. Here, r r r = r r r + r r r . The coefficients of the vector addition theorem A νm ,nm and B νm ,nm can be expressed in terms of summation over Wigner 3j symbols [START_REF] Abramowitz | Handbook of Mathematical Functions: With Formulas, Graphs, and Mathematical Tables[END_REF], calculating which involve computations of a large number of factorials. For the configuration shown in Fig. 1 (with m = m due to axisymmetry condition), the vector translation coefficients are given by (using k = k f for the coordinate systems O and O positioned in free-space):

A νm,nm = ν(ν + 1) n(n + 1) 2π ν(ν + 1) ν+n ∑ p=|ν-n| n(n + 1) + ν(ν + 1) -p(p + 1) × i ν+p-n ψ p,0 (k f r ) A(m, n, -m, ν, p) , (11)
B νm,n,m = ν(ν + 1) n(n + 1) imk f r ν(ν + 1) ν+n ∑ p=|ν-n| i ν+p-n ψ p,0 (k f r ) A(m, n, -m, ν, p), (12)
where,

A(m, n, -m, ν, p)=(-1) m (2n+1)(2ν +1)(2p+1) 4π × n ν p 0 0 0 n ν p -m m 0 , (13)
with the symbol n ν p m µ q representing the Wigner 3j cooefficient. The above relations differ from those that are given in Ref. [START_REF] Chew | Derivation of the vector addition theorem[END_REF] due to the different definitions of M M M (p) nm (kr r r) and N N N (p) nm (kr r r) waves that we are using (given in Eqs. ( 6) and ( 7)). The corresponding scalar coefficient β νm,nm is given by (Eq. 25 in Ref. [START_REF] Chew | Derivation of the vector addition theorem[END_REF] modified to reflect axisymmetry condition for our configuration):

β νm,nm (k f r ) = ν+n ∑ p=|ν-n| 4π i ν+p-n z (3) p (k f r ) 2p + 1 4π A(m, n, -m, ν, p). (14)
In Eqs. ( 11), [START_REF] Loomis | Theory of heat transfer by evanescent electromagnetic waves[END_REF], and ( 14), the summation over p occurs in steps of 2, i.e., p

= |ν -n|, |ν - n| + 2, • • • , ν + n -2, ν + n.
When large number of coefficients are required (an expression for number of terms required for convergence in the summations of Eqs. ( 9) and ( 10) for radiative energy transfer calculations between two spheres has been derived in Ref. [START_REF] Sasihithlu | Convergence of vector spherical wave expansion method applied to nearfield radiative transfer[END_REF]. We will show in Sec. 4 that n and ν ∼ 10 3 when d/R 2 ∼ 10 -3 ) computation of these Wigner 3j symbols becomes computationally tedious. Recursion relations derived for the scalar addition theorem [START_REF] Chew | Recurrence relations for three-dimensional scalar addition theorem[END_REF] and the vector addition theorem [START_REF] Chew | Efficient ways to compute the vector addition theorem[END_REF] can be used, resulting in reduced computational times. The recursion relation for the scalar translation coefficient β is given by (Eq. 23 and Eq. 27 in Ref. [START_REF] Chew | Recurrence relations for three-dimensional scalar addition theorem[END_REF]):

β νm,n+1,m = 1 a + nm -a - nm β νm,n-1,m + a + ν-1,m β ν-1 m,nm + a - ν+1,m β ν+1 m,nm , (15)
and

β ν n+1,n+1 n+1 = 1 b + nn b + ν-1,n β ν-1 n,nn + b - ν+1,n β ν+1 n,nn . (16)
The form of the constants a + nm , a - nm and b + nm , b - nm have been retained from Ref. [START_REF] Chew | Recurrence relations for three-dimensional scalar addition theorem[END_REF]. Since the vector translation coefficients A νm,nm and B νm,nm can be related to β , they are computed directly from β as (Eq. 12a and Eq. 12b in Ref. [START_REF] Chew | Efficient ways to compute the vector addition theorem[END_REF] changed according to the axisymmetry of our configuration and definition of M M M and N N N waves):

A νm,nm = ν(ν + 1) n(n + 1) β νm,nm + k f r ν + 1 (ν + m + 1)(ν -m + 1) (2ν + 1)(2ν + 3) β ν+1,m,nm + k f r ν (ν + m)(ν -m) (2ν -1)(2ν + 1) β ν-1,m,nm , (17)
and

B νm,n,m = imk f r ν(ν + 1)n(n + 1) β νm,n,m . (18)
To initialize the recursion for the scalar translation coefficient the following expression is used :

β ν0,00 = (-1) ν √ 2ν + 1 z (3) ν (k f r ), (19)
where z

ν (k f r ) represents the spherical Hankel function. When ν k f r , spherical Hankel and Bessel functions can be approximated by the following asymptotic forms:

z (3) ν (k f r ) ≈ i 4 (2ν + 1)k f r 2ν + 1 ek f r (ν+1/2) , (20) z (1)
ν (k f r ) ≈ 4 (2ν + 1)k f r ek f r 2ν + 1 (ν+1/2) . ( 21
)
It follows from Eq. ( 20) and Eq. ( 21) that evaluating spherical Hankel and Bessel functions in the limit ν k f r can lead to numbers too large or too small for a given floating point format on a computer. For instance, for typical values of ν = 1000 and k f r = 200 that we encounter, z

ν (k f r ) ≈ 8.3 × 10 567 and z

(1) ν (k f r ) ≈ 3.05 × 10 -574 , while the maximum and minimum positive floating point numbers for double precision floating point format are 1.7976 × 10 +308 and 2.225 × 10 -308 respectively (both approximate). However, since products of such large and small terms are important for subsequent calculations (see for example, Eq. ( 25)), we need to modify the algorithm so that evaluations such large or small numbers can be minimized. This can be achieved through the "normalization" of the translation coefficients with appropriate factors, which we explain below.

A careful look at the form of the spectral radiative conductance (Eq. 28 in Ref. [START_REF] Narayanaswamy | Thermal near-field radiative transfer between two spheres[END_REF]) gives us a hint about the appropriate factors. It is observed that coefficients D lM νm and D lN νm which represent the coefficients of M M M and N N N waves in the scattered field in vacuum (Eq. 19 in Ref. [START_REF] Narayanaswamy | Thermal near-field radiative transfer between two spheres[END_REF]) are accompanied by factors of the form z

(1) l (k a R 1 )/z (1) β (k f R 2 ) and z (1) l (k a R 1 )/ζ (1) β (k f R 2 ) respec- tively.
Hence, what is important for the computation of the radiative energy transfer are not the terms D lM νm and D lN νm themselves but D lM νm z

(1)

l (k a R 1 )/z (1) ν (k f R 2 ) and D lN νm z (1)
l (k a R 1 )/ζ (1)
ν (k f R 2 ). These factors can be used in the coupled linear equations for the coefficients of vector spherical waves in vacuum (Eq. 20 in Ref. [START_REF] Narayanaswamy | Thermal near-field radiative transfer between two spheres[END_REF]) to arrive at different forms of the normalizing factors for the translation coefficients A νm,nm and B νm,nm as has been shown in Appendix 5.

Consider one representative form for our discussion: z

ν (k f R 1 )A νm,nm /z (1)
n (k f R 2 ). Since A νm,nm and B νm,nm are related to β νm,nm as given in Eqs. ( 17) and [START_REF] Domingues | Heat transfer between two nanoparticles through near-field interaction[END_REF], it is sufficient if we build recursion relations for the normalized scalar translation coefficient β N νm,nm given by Figures 3(a), 3(b), 3(c) and 3(d) show the contour plots of of logarithm of the expression in Eq. ( 31) for two spheres of radius ratio R 2 /R 1 = 1, 3, 10, and 20, respectively, as a function of n and ν. The radius of the smaller sphere, R 1 , and the minimum gap between the spheres, d, are taken to be 10 µm and 50 nm for all the cases. When R 2 /R 1 = 1 it can be observed from Fig. 3(a) that the dominant terms in the matrix are present in a band along the diagonal with terms on either side of the band decreasing exponentially. This behavior has been utilized to simplify the computation of translation coefficients and the linear equations arising from imposing boundary conditions. Terms in the matrix with absolute value less than 10 -6 can be approximated to be zero and sparse routines can be employed to solve the linear equations. This approximation results in an error of less than 0.05 % in the final computed value of the conductance. For R 2 /R 1 > 1 the distribution of the dominant terms in the matrix changes markedly. From Fig. 3(c) and 3(d), where R 2 /R 1 1, it can be observed that only the first few rows of the matrix need be computed since the absolute values of the terms in the higher rows fall off exponentially. This substantially frees up memory resources while also reducing the computation time considerably. The errors from the two approximations (one term approximation, and making use of only the dominant (n, ν) terms in the normalized translation coefficient β N νm,nm ) are compared in Table 1.

The variation of β νm,nm with R 2 /R 1 can be determined from the asymptotic behavior of the expression in Eq. [START_REF] Narayanaswamy | Near-field radiative heat transfer between a sphere and a substrate[END_REF], where

k f r = k f R 1 + k f R 2 + δ with δ = k f d → 0 and n k f R 1 and ν k f R 2 .
From Ref. [START_REF] Abramowitz | Handbook of Mathematical Functions: With Formulas, Graphs, and Mathematical Tables[END_REF], the asymptotic form for z

(1)

ν (k f R 1 ) for ν k f R 1 is given by : z (1) ν (k f R 1 ) ≈ exp (ν +1/2) 2 -(k f R 1 ) 2 -(ν +1/2) cosh -1 [ (ν+1/2) k f R 1 ] 2 k f R 1 (ν +1/2) 2 -(k f R 1 ) 2 . ( 32
)
Similar asymptotic forms for z

n (k f R 2 ) and z

ν (k f R 1 ) z (3) n+ν (k f r )/z (3) n (k f R 2 )
| as a function of n and ν for two spheres with successive radius ratios R 2 /R 1 = (a) 1, (b) 3, (c) 10, and (d) 20 with R 1 = 10 µm, and the minimum gap maintained at 50 nm for all the cases. The dashed-lines denotes the contour line for a value of -6 which is taken as the cutoff point below which values for the normalized vector translation coefficients are approximated to zero. The line of maximum (shown as dotted lines) given by Eq. ( 35) has been superimposed on these contour plots or n+ν. When n 1/2, ν 1/2, the magnitude of the expression in Eq. ( 31) is approximately:

exp ν 2 -(k f R 1 ) 2 + n 2 -(k f R 2 ) 2 -(n + ν) 2 -(k f R 1 + k f R 2 + δ ) 2 × exp (n + ν) cosh -1 n+ν k f R 1 +k f R 2 +δ -n cosh -1 n k f R 2 -ν cosh -1 ν k f R 1 . ( 33
)
When δ reduces to zero, employing cosh -1 x = log(x + √ x 2 -1), taking n k f R 1 and ν is the spectral conductance value without any approximations, G A 1 ω is the spectral conductance value obtained by using the one-term approximation and G A 2 ω is obtained by considering the dependence of normalized translation coefficients on the radius ratio of the two spheres. All conductances are in units of nW.K -1 . The errors from these two approximations are denoted by 'Error A1' and 'Error A2' respectively (both in %). The values shown here are computed at the surface phonon-polariton frequency 0.061 eV for two cases (1) R 2 = R 1 , and (2) R 2 = 10R 1 . R 1 is kept a constant at 13.7 µm. k f R 2 and rearranging, the expression in Eq. ( 33) can be written as:

Case

G Exact ω G A 1 ω Error A1 G A 2 ω Error A2 R 2 =
1 + n ν n 1 + ν n ν 1 + R 2 R 1 -ν 1 + R 1 R 2 -n . ( 34
)
The expression in Eq. ( 34) can be shown to attain a maximum value of 1 when

n ν = R 2 R 1 ( 35
)
This condition conforms well with what is observed in the contour plots of Fig. 3 where the line representing Eq. ( 35) has been superimposed on each of the contour plots. Further analysis of the expression in Eq. ( 31) provides us with an estimate of the number of terms that are important along line of maximum. Using the following Taylor series expansion about δ = 0 (truncating it to the first order) and substituting Eq. ( 35), we can see that the expression in Eq. ( 31) simplifies to:

exp - δ ν 2 -(k f R 1 ) 2 ν 2 k f R 1 -k f R 1 . ( 36
)
The above expression becomes negligibly small when δ ν

2 -(k f R 1 ) 2 k f R 1 , i.e., if ν R 1 d 2 + (k f R 1 ) 2 , (37)
From Eq. ( 35) the condition for n would be :

n R 2 R 1 R 1 d 2 + (k f R 1 ) 2 . ( 38
)
This condition is similar to the criterion for the convergence of the vector spherical wave expansion method of the near-field radiative transfer between two spheres which was shown in Ref. [START_REF] Sasihithlu | Convergence of vector spherical wave expansion method applied to nearfield radiative transfer[END_REF] to be of the form:

N conv = 2 R 2 d + ek f r 2 , (39)
This criterion was derived in Ref. [START_REF] Sasihithlu | Convergence of vector spherical wave expansion method applied to nearfield radiative transfer[END_REF] heuristically based on comparison with the convergence criterion for the near-field radiative transfer between planar bodies. Equations ( 37) and ( 38) 41)) to predict the far-field contribution to the conductance demonstrate explicitly the role of translation coefficients in attaining convergence of near-field radiative transfer calculations between two spheres. For our numerical calculations we have used Eq. ( 39) as the criterion for the number of vector spherical waves since it is typically of higher value than that represented in Eq. ( 38).

Validity of modified proximity approximation

The conductance G can be calculated once the DGF for the two sphere configuration is determined using the theory developed in Ref. [START_REF] Narayanaswamy | Thermal near-field radiative transfer between two spheres[END_REF], modified suitably to account for the changes described in Secs. 2-4. The radiative conductance G (units WK -1 ) between the two spheres shown in Fig. 1 is defined as:

G = lim T A →T B P (T A , T B ) |T A -T B | , (40)
where T A and T B are the temperatures of the two spheres, and P (T A , T B ) is the rate of heat transfer between them. As explained in the introduction, it is also possible to predict the radiative transfer between nonplanar surfaces by using the MPA, i.e., by treating the nonplanar surfaces as a series of parallel surfaces with varying gaps, as shown in Fig. 4(a). The effectiveness of MPA in predicting near-field radiative conductance for two equally sized spheres has been shown in [START_REF] Sasihithlu | Proximity effects in radiative heat transfer[END_REF]. Here, we investigate the relevance of MPA for the case when R 2 R 1 .

Using the MPA, the conductance G MPA between two unequal sized spheres of radii R 1 and R 2 for a minimum separation gap d and temperature T can be written as:

G MPA (d, T ) = R 1 0 h n f (z)2πr dr + G c (d, T ), (41)
where R 1 is the radius of the smaller sphere,

z = d + R 1 + R 2 -R 2 1 -r 2 -R 2 2 -r 2
is the gap at distance r from the symmetry axis as shown in Fig. 4(a), h n f (z) is the near-field contribution [START_REF] Sasihithlu | Proximity effects in radiative heat transfer[END_REF] to the radiative heat transfer coefficient between two half-spaces at gap z and G c (d, T ) can be approximated by the conductance value from classical radiative transfer theory when diffraction effects are negligible. G c (d, T ) for two unequal spheres of equal emissivity ε and temperature T is given by [START_REF] Siegel | Thermal Radiation Heat Transfer[END_REF]:

G c (d, T ) = 4σ T 3 (4πR 2 1 ) [(1 -ε)/ε](1 + R 2 1 /R 2 2 ) + 1/F 12 (d) , (42)
where F 12 is the gap dependent view factor between the two spheres and can be approximated (to an accuracy of 1%) by [START_REF] Juul | Investigation of approximate methods for calculation of the diffuse radiation configuration view factor between two spheres[END_REF]:

F 12 (d) = 1 2 1 -1 - 1 (d/R 2 + R 1 /R 2 + 1) 2 . ( 43
)
The emissivity for a silica half-space has been computed and plotted as a function of frequency in eV in Fig. 4(b). The conductance values have been computed by integrating over the frequency range 0.041 eV to 0.164 eV, and the value of G c (d, T ) has been appropriately adjusted to reflect this [43, see supplemental information]. For the radius of the spheres that we have considered in our study (R 1 = 2.5 µm, 13.7 µm) the value of G(d) at d ≈ 40 µm (when near-field effects are negligible) is taken to be the value G c at that gap and the effect of gap-dependence of view factor for smaller gaps is included by using:

G c (d 1 , T ) G c (d 2 , T ) = [(1 -ε)/ε](1 + R 2 1 /R 2 2 ) + 1/F 12 (d 2 ) [(1 -ε)/ε](1 + R 2 1 /R 2 2 ) + 1/F 12 (d 1 ) , (44)
The comparison between G and G MPA is shown in Fig. 5(a) for two spheres with R 1 = 13.7 µm and 2.5 µm and R 2 = 40R 1 . The error between G and G MPA is plotted in Fig. 5(b). For gaps d/R 1 < 0.1, MPA is observed to be able to predict the exact computed values of the conductance with errors less than 1%.

Since the total radiative conductance between the two spheres has contributions from frequencies where there is resonant enhancement from surface phonon polaritons as well as contributions from non-resonant frequencies where surface phonon polaritons are not active (see Fig. 2 in [START_REF] Sasihithlu | Proximity effects in radiative heat transfer[END_REF]), it would be of interest to observe if the MPA accurately predicts the contribution from both these regions. Computed values (by DGF formalism) and MPA predictions of the spectral conductance G ω at a resonant (0.061 eV) and a non-resonant (0.081 eV) frequency are plotted in Fig. 6(a) and Fig. 7(a). The contribution from the far-field radiative conductance G c (d, T ) in Eq. [START_REF] Chew | Derivation of the vector addition theorem[END_REF] has to be appropriately modified to reflect spectral conductance. From proximity approximation theory, as d/R 1 → 0, the spectral conductance G ω (d) for two spheres of unequal radii R 1 and R 2 is expected to scale as [START_REF] Sernelius | Surface Modes in Physics[END_REF] :

G ω (d) ∼ R 1 R 2 R 1 + R 2 1 d ≈ R 1 d (for R 2 R 1 ). (45)
This characteristic R 1 /d behavior is observed for the resonant frequency contributions shown in Fig. 6(a) for d/R 1 0.02 where a slope of -1 is indicated. However such behavior is not observed for the non-resonant frequency contributions shown in Fig. 7(a). The error between G ω and the spectral conductance predicted by MPA, G MPA ω , for the resonant and non-resonant frequency contributions are shown in Fig. 6(b) and Fig. 7(b) respectively. In the far-field region (d/R 1 2 for R 1 = 13.7 µm, and d/R 1 10 for R 1 = 2.5 µm) where the enhancement due to tunneling of waves (surface waves at the resonant frequency and evanescent waves at non-resonant frequency) is negligible, the form of MPA in Eq. ( 41) predicts that the variation in G ω with gap is primarily due to the changing view factor between the spheres d/R 1 0.07, for R 1 = 13.7 µm) the variation of G ω with gap is dependent on both the changing view factor with gap as well as increased tunneling of waves. For gaps d/R 1 0.07 (below which the view factor increases by less than 1%) the enhanced radiative transfer with decreasing gap can be attributed entirely due to increased tunneling of waves. At such small gaps MPA is able to model the enhancement at the resonant frequency within ≈ 5% errors irrespective of the value of R 1 , whereas at the non-resonant frequency the error is observed to be greater than 10% when R 1 = 2.5 µm. Despite such high errors at the non-resonant frequencies, MPA is successful in predicting the overall conductance when R 1 = 2.5 µm with error less than 1% for d/R 1 0.1 as observed in Fig. 5. This apparent discrepancy has been analyzed in detail in Ref. [START_REF] Gu | Proximity approximation and radiative transfer between sub-wavelength spheres and planar surfaces[END_REF] and has been attributed to decreased contribution from non-resonant frequencies as R 1 decreases.

We summarize the main contributions from this work: (a) The exponential behavior of the spherical Hankel function z

n (k f r ) when n ≥ k f r has been utilized to propose a simplified form, referred to in this work as the one-term approximation, for the coefficients in the translation addition theorem. The one-term approximation is valid when n ≥ Ck f r where C is a constant which is chosen depending on the desired accuracy. The recursion relations for these simplified translation forms are also given. They are simpler (b) A method to normalize the translation coefficients has been proposed and the behavior of the normalized translation coefficients on the relative size of the two spheres has been exploited to simplify the computations of near-field radiative transfer between two spheres.

(c) These simplifications are utilized to compute near-field radiative transfer between a mesoscopic and a macroscopic sphere with size ratio of 40 and we show that there is agreement with error less than 1% between the predictions of the modified proximity approximation and the exact computations for d/R 1 < 0.1. Thus the modified proximity approximation, which was initially proposed for the near-field radiative transfer between two equal sized spheres, has been shown to be valid for a more general configuration of unequal sized spheres as well.

( 1 )

 1 l (kr) is the spherical Bessel function of order l. For j = 3, the M M M and N N N waves are outgoing spherical waves that are singular at the origin and z

( 3 )

 3 n (kr) is the spherical Hankel function of the first kind of order n. The radial function ζ

(

Fig. 1 .

 1 Fig. 1. The configuration for this study consisting of two spheres of unequal radii R 1 and R 2 (labeled sphere 'a' and sphere 'b' respectively) and the relation between r r r , r r r and r r r . The surface to surface gap between the two spheres d is given by d = r -R 1 -R 2

( 3 )Fig. 3 .

 33 Fig. 3. Contour plots of the expression log 10 |z (1)

Fig. 4 .

 4 Fig. 4. (a) Application of MPA for the two spheres where the curved are approximated by a series of flat surfaces with varying gaps z. (b) the plot of spectral emissivity for a silica half-plane as a function of frequency in eV which is used in the form of MPA (Eq. (41)) to predict the far-field contribution to the conductance

Fig. 5 .

 5 Fig. 5. (a) Plot of computed values of the total conductance (dotted) and the MPA (solid line)as a function of the non-dimensional gap d/R 1 for two spheres with R 2 /R 1 = 40 . The study has been performed for R 1 = 13.7µm and 2.5µm (b) The % error between the computed values and values from the MPA are plotted as a function of d/R 1

Fig. 6 .

 6 Fig. 6. (a) Comparison between the computed values of G ω at a surface phonon-polariton frequency (0.061 eV) (dotted) and G MPA ω (solid line) as a function of the non-dimensional gap d/R 1 for two spheres of radius R 1 = 2.5 µm, 13.7 µm and R 2 = 40R 1 (b) The % error between G ω and G MPA ω as a function of d/R 1

Fig. 7 .

 7 Fig. 7. (a) Comparison between the computed values of G ω at frequency 0.0801 eV (dotted) and G MPA ω (solid line) as a function of the non-dimensional gap d/R 1 for two spheres of radius R 1 = 2.5µm, 13.7µm and R 2 = 40R 1 (b) The % error between G ω and G MPA ω as a function of d/R 1

Table 1 .

 1 Table showing the error from the two approximations. Here, G Exact ω

Acknowledgment

This work is funded partially by ONR Grant N00014-12-1-0996.

From Eq. ( 22) and Eq. ( 15), Eq. ( 16) it can be shown that the recursion relations to compute β N νm,nm are given by :

n-1 (k f R 1 )

and

These recursion relations are only dependent on ratio of successive orders of spherical Bessel/Hankel functions and, hence, we do not encounter very small or large numbers when ν k f r . The ratios have been computed by method of continued fractions [START_REF] Cuyt | Handbook of Continued Fractions for Special Functions[END_REF]. The initialization for this recursion corresponding to Eq. ( 19) will be:

The spherical Hankel functions in Eq. ( 25), if computed individually, will lead to diverging numbers when ν k f r as explained previously. To circumvent this, β N ν0,00 in Eq. ( 25) is determined by computing the logarithm of each of the spherical Bessel and Hankel functions. The procedure to compute logarithm of Bessel functions and Hankel functions is detailed in Ref. [START_REF] Rothwell | Computation of the logarithm of bessel functions of complex argument[END_REF].

In addition to circumventing the possibility of encountering diverging numbers, the usage of normalized translation coefficients enables further simplification of the computation procedure as explained in Sec. 4.

One term approximation

In most cases arrays of A νm,nm and B νm,nm are needed for n = 1, ..., n max and ν = 1, ..., ν max at fixed values of m. As mentioned in Ref. [START_REF] Chew | Efficient ways to compute the vector addition theorem[END_REF] when n max = ν max = 10, the speed when using the recurrence relations given in Eq. ( 15) is about 400 times faster than using the Wigner 3j coefficients in Eq. [START_REF] Mulet | Enhanced radiative transfer at nanometric distances[END_REF]. However, when computing the near-field radiative transfer between two spheres separated by a small gap, n max and ν max can easily be as high as 1000 and even the recurrence relations mentioned above become computationally expensive. We can show that when the number of spherical waves needed for convergence n max , ν max k f r , the recurrence relations can be simplified further to make them more computationally viable.

Consider the form for the scalar translation coefficient β νm,nm (k f r ) given in Eq. ( 14). We note that the expression for β νm,nm (k f r ) depends on a summation over different orders p of the spherical Hankel function z 20)) it is possible to show that for p k f r the contribution from the last term in the summation i.e p = n + ν dominates over the rest of the terms. Hence, only the term with p = n + ν need be retained in the summation. From Eq. ( 14), the one-term approximation for the scalar translation coefficient β νm,nm can be written as :

An estimate of the error from retaining only the last term in Eq. ( 14) can be obtained by comparing the magnitudes of the p = n + ν and p = n + ν -2 terms. Since A(m, n, -m, ν, p) and √ 2p + 1 vary relatively weakly with p, in comparison to z

p (k f r ), the error from retaining only the term p = n + ν in the summation in Eq. ( 14) can be estimated to be:

Using Eq. ( 20) in the limit p 1 we can show that :

The error from Eq. ( 27) then reduces to

For 1 % error we get n + ν ≈ 5k f r . For our computations we have used n + ν ≥ 7k f r as the criterion for employing the one-term approximation. From Fig. 2 it can be observed that when

n+ν-2 (k f r ). The error (in %) in the spectral conductance at the surface phonon-polariton frequency, which is ≈ 0.061 eV for silica, (see Fig. 2 in Ref. [START_REF] Sasihithlu | Proximity effects in radiative heat transfer[END_REF]) when the one-term approximation for β νm,nm (k f r ) is used beyond n + ν = pk f r is shown in Fig. 2(b) for different values of n + ν. The errors are computed with respect to the spectral conductance value when the approximated form for translation coefficients is not employed. It can be observed that adopting n + ν = 7k f r as the criterion for employing the one-term approximation gives us an error of about 0.05 % in the spectral conductance.

n+ν-2 (k f r )) as a function of n+ν for k f r = 33 (arbitrarily chosen). The point n + ν = 7 k f r , beyond which the one-term approximation has been adopted in our computations for calculating the translation coefficients, is marked in in the figure . (b) The error in spectral conductance at the surface phonon-polariton frequency (0.061 eV) when different values of (n + ν)/k f r are chosen as the criteria for employing the one-term approximation. From the plot, a criterion n + ν = 7k f r is observed to give an error of ≈ 0.02% in the spectral conductance. The spectral conductance has been computed for two spheres of size R 1 = 13.7 µm and R 2 = 40R 1 with minimum gap d/R 1 = 0.01.

Since computation of Wigner 3j coefficient is computationally tedious, recursion relations can be derived for the approximated form for β νm,nm . Employing

it is easy to show that the recursion relations turn out to be :

The form of recursion given in Eq. ( 30) is computationally simpler than the corresponding recursion relation for the exact form of β νm,nm given in Eq. [START_REF] Pendry | Radiative exchange of heat between nanostructures[END_REF]. Once the scalar translation coefficients β νm,nm are computed, the vector translation coefficients A νm,nm and B νm,nm can be got from Eqs. ( 17) and ( 18).

Dependence of normalized translation coefficients on the radius ratio of the spheres

For near-field scattering problems when the closest gap between the two spheres d R 1 , R 2 (here, d = r -R 1 -R 2 ), the normalized translation coefficients β N nm,νm show a marked dependence on the radius ratio R 2 /R 1 . This behavior can be used to improve the efficiency of computation of near-field quantities considerably. In the one-term approximation, β N νm,nm varies with n and ν mainly through the ratio:

Appendix: Normalizing factors for the vector translation coefficients A νm,nm and B νm,nm

The coupled linear equations are given by (Eq. 20 in Ref. [START_REF] Narayanaswamy | Thermal near-field radiative transfer between two spheres[END_REF]):

where the symbol (m, 1) denotes the maximum of m and 1, u n (R 1 ) and v n (R 1 ) are given by:

and

ν (k f R 2 ) respectively. The equivalent normalizing factors for C lM nm and C lN nm are z

(1)

n (k f R 1 ) respectively. Using these in Eqs. ( 46)-( 49) and rearranging, the coupled linear equations reduce to:

A νm,nm (-k f D)

l (k a R 1 ) z

(1)

B νm,nm (+k f D)

A νm,nm (+k f D) ν (k f R 1 ) for ν k f R 1 only one of them (the latter) has been chosen as a representative form for discussion in Eq. [START_REF] Krüger | Nonequilibrium electromagnetic fluctuations: Heat transfer and interactions[END_REF].