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

Karthik Sasihithlu; Arvind Narayanaswamy

doi:10.1364/OE.22.014473

1. 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 [1–4], thermal rectification [5], nanopatterning [6], thermally assisted magnetic recording [7], thermal imaging [8] and noncontact radiative cooling [9]. The theoretical analysis has relied on Rytov’s theory of fluctuational electrodynamics [10] and has been used to compute near-field radiative transfer between two flat surfaces [11–16], between two dipoles [17, 18], and between a dipole and a flat surface [15, 19]. Only recently has this been applied to other more experimentally feasible geometries such as between two mesoscopic spheres [20] and between a mesoscopic sphere and a planar surface [21, 22]. 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 ≈ ħc/k_BT, where 2πħ 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 [22–26]. Recent experimental validation of this near-field enhancement for sub-micron gaps has been limited to that between a mesoscopic sphere and a plate [27–31] and between a scanning probe tip and a planar surface [32–34].

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) [35, 36] 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 [37–39]. This procedure has been adopted to compute the thermal near-field radiative transfer between two spheres [20] and between a sphere and a plane [21].

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 [40, 41], recursion relations have to be used to ensure their efficient computation [42].

In addition to the exact method [20] 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) [43] differs from the proximity force approximation (PFA) [44, 45] 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 [46] and the near-field radiative transfer [47] 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 [43]. Hence it would be of interest to also study the relevance of the MPA when the radius of one of the spheres R₂ is increased such that R₂ ≫ R₁ where R₁ is the radius of the smaller sphere, since when R₂/R₁ → ∞ 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 [20], 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.

2. 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 z-axis as shown in Fig. 1. Emission and scattering of electromagnetic waves from such a system involves solving the vector Helmholtz equation given by:

\nabla \times \nabla \times P (r) - k^{2} P (r) = 0,

where P(r) is the electric or magnetic field at position vector 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 [20]:

M_{n m}^{(j)} (k r) = z_{n}^{(j)} (k r) V_{n m}^{(2)} (θ, ϕ),

N_{n m}^{(j)} (k r) = ζ_{n}^{j} (k r) V_{n m}^{(3)} (θ, ϕ) + \frac{z_{n}^{(j)} (k r)}{k r} \sqrt{n (n + 1)} V_{n m}^{(1)} (θ, ϕ),

where

M_{n m}^{(j)} (k r)

and

N_{n m}^{(j)} (k 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 and N waves are regular waves and remain finite at the origin and

z_{l}^{(1)} (k r)

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

z_{n}^{(3)} (k r)

is the spherical Hankel function of the first kind of order n. The radial function

ζ_{n}^{j} (x) = \frac{1}{x} \frac{d}{d x} (x z_{n}^{(j)} (x))

.

V_{n m}^{(1)} (θ, ϕ)

,

V_{n m}^{(2)} (θ, ϕ)

, and

V_{n m}^{(3)} (θ, ϕ)

are vector spherical harmonics of order (n, m) and are given by:

V_{n m}^{(1)} (θ, ϕ) = \hat{r} Y_{n m},

V_{n m}^{(2)} (θ, ϕ) = \frac{1}{\sqrt{n (n + 1)}} (- \hat{ϕ} \frac{\partial Y_{n m}}{\partial θ} + \hat{θ} \frac{i m}{\sin θ} Y_{n m}),

and

V_{n m}^{(3)} (θ, ϕ) = \frac{1}{\sqrt{n (n + 1)}} (\hat{θ} \frac{\partial Y_{n m}}{\partial θ} + \hat{ϕ} \frac{i m}{\sin θ} Y_{n m}) .

Fig. 1 The configuration for this study consisting of two spheres of unequal radii R₁ and R₂ (labeled sphere ‘a’ and sphere ‘b’ respectively) and the relation between r, r′ and r″. The surface to surface gap between the two spheres d is given by d = r″ − R₁ − R₂

Download Full Size | PDF

Following Ref. [42] 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 ∇²ψ + k²ψ = 0 are given by:

ψ_{n m}^{(j)} (k, r) = z_{n}^{(j)} (k r) Y_{n m} (θ, ϕ),

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_{n m}^{(j)} (k, r) = \frac{1}{\sqrt{n (n + 1)}} \nabla \times r ψ_{n m}^{(j)} (k, r),

and

N_{n m}^{(j)} (k, r) = \frac{1}{\sqrt{n (n + 1)}} \frac{1}{k} \nabla \times \nabla \times r ψ_{n m}^{(j)} (k, r) .

The representations for M and N waves in Eqs. (2), (3), (6), and (7) are equivalent. They differ from the corresponding representations in Ref. [42] in only a factor

1 / \sqrt{n (n + 1)}

which has been introduced to ensure that

V_{n m}^{(1)} (θ, ϕ)

,

V_{n m}^{(2)} (θ, ϕ)

, and

V_{n m}^{(3)} (θ, ϕ)

are orthonormal to each other.

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

ψ_{n m} (k, r) = \sum_{ν = 1}^{\infty} \sum_{m^{'} = - ν}^{ν} ψ_{ν m^{'}} (k, r^{'}) β_{ν m^{'}, n m},

and the vector addition theorem for the vector spherical waves

M_{n m}^{(j)}

and

N_{n m}^{(j)}

waves is given by:

M_{n m}^{(j)} (k r) = \sum_{ν = 1}^{\infty} \sum_{m^{'} = - ν}^{ν} [M_{ν m^{'}}^{(j)} (k r^{'}) A_{ν m^{'}, n m} + N_{ν m^{'}}^{(j)} (k r^{'}) B_{ν m^{'}, n m}],

N_{n m}^{(j)} (k r) = \sum_{ν = 1}^{ν = \infty} \sum_{m^{'} = - ν}^{ν} [N_{ν m^{'}}^{(j)} (k r^{'}) A_{ν m^{'}, n m} + M_{ν m^{'}}^{(j)} (k r^{'}) B_{ν m^{'}, n m}]

where r and r′ are coordinates with respect to the coordinate systems centered at O and O′ respectively as shown in Fig. 1. Here, 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 [48], 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):

\begin{matrix} A_{ν m, n m} = \sqrt{\frac{ν (ν + 1)}{n (n + 1)}} {\frac{2 π}{ν (ν + 1)} \sum_{p = | ν - n |}^{ν + n} [n (n + 1) + ν (ν + 1) - p (p + 1)] \times \\ i^{ν + p - n} ψ_{p, 0} (k_{f} r^{″}) A (m, n, - m, ν, p)}, \end{matrix}

B_{ν m, n, m} = \sqrt{\frac{ν (ν + 1)}{n (n + 1)}} \frac{i m k_{f} r^{″}}{ν (ν + 1)} \sum_{p = | ν - n |}^{ν + n} i^{ν + p - n} ψ_{p, 0} (k_{f} r^{″}) A (m, n, - m, ν, p),

where,

A (m, n, - m, ν, p) = {(- 1)}^{m} \sqrt{\frac{(2 n + 1) (2 ν + 1) (2 p + 1)}{4 π}} \times (\begin{matrix} n & ν & p \\ 0 & 0 & 0 \end{matrix}) (\begin{matrix} n & ν & p \\ - m & m & 0 \end{matrix}),

with the symbol

(\begin{matrix} n & ν & p \\ m & μ & q \end{matrix})

representing the Wigner 3j cooefficient. The above relations differ from those that are given in Ref. [41] due to the different definitions of

M_{n m}^{(p)} (k r)

and

N_{n m}^{(p)} (k 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. [41] modified to reflect axisymmetry condition for our configuration):

β_{ν m, n m} (k_{f} r^{″}) = \sum_{p = | ν - n |}^{ν + n} 4 π i^{ν + p - n} z_{p}^{(3)} (k_{f} r^{″}) \sqrt{\frac{2 p + 1}{4 π}} A (m, n, - m, ν, p) .

In Eqs. (11), (12), 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. [49]. We will show in Sec. 4 that n and ν ∼ 10³ when d/R₂ ∼ 10⁻³) computation of these Wigner 3j symbols becomes computationally tedious. Recursion relations derived for the scalar addition theorem [50] and the vector addition theorem [42] 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. [50]):

β_{ν m, n + 1, m} = \frac{1}{a_{n m}^{+}} (- a_{n m}^{-} β_{ν m, n - 1, m} + a_{ν - 1, m}^{+} β_{ν - 1 m, n m} + a_{ν + 1, m}^{-} β_{ν + 1 m, n m}),

and

β_{ν n + 1, n + 1 n + 1} = \frac{1}{b_{n n}^{+}} (b_{ν - 1, n}^{+} β_{ν - 1 n, n n} + b_{ν + 1, n}^{-} β_{ν + 1 n, n n}) .

The form of the constants

a_{n m}^{+}

,

a_{n m}^{-}

and

b_{n m}^{+}

,

b_{n m}^{-}

have been retained from Ref. [50].

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. [42] changed according to the axisymmetry of our configuration and definition of M and N waves):

\begin{matrix} A_{ν m, n m} = \sqrt{\frac{ν (ν + 1)}{n (n + 1)}} {β_{ν m, n m} + \frac{k_{f} r^{″}}{(ν + 1)} \sqrt{\frac{(ν + m + 1) (ν - m + 1)}{(2 ν + 1) (2 ν + 3)}} β_{ν + 1, m, n m} + \\ \frac{k_{f} r^{″}}{ν} \sqrt{\frac{(ν + m) (ν - m)}{(2 ν - 1) (2 ν + 1)}} β_{ν - 1, m, n m}}, \end{matrix}

and

B_{ν m, n, m} \frac{i m k_{f} r^{″}}{\sqrt{ν (ν + 1) n (n + 1)}} β_{ν m, n, m} .

To initialize the recursion for the scalar translation coefficient the following expression is used :

β_{ν 0, 00} = {(- 1)}^{ν} \sqrt{2 ν + 1} z_{ν}^{(3)} (k_{f} r^{″}),

where

z_{ν}^{(3)} (k_{f} r^{″})

represents the spherical Hankel function. When ν ≫ k_fr″, spherical Hankel and Bessel functions can be approximated by the following asymptotic forms:

z_{ν}^{(3)} (k_{f} r^{″}) \approx i \sqrt{\frac{4}{(2 ν + 1) k_{f} r^{″}}} {(\frac{2 ν + 1}{e k_{f} r^{″}})}^{(ν + 1 / 2)},

z_{ν}^{(1)} (k_{f} r^{″}) \approx \sqrt{\frac{4}{(2 ν + 1) k_{f} r^{″}}} {(\frac{e k_{f} r^{″}}{2 ν + 1})}^{(ν + 1 / 2)} .

It follows from Eq. (20) and Eq. (21) that evaluating spherical Hankel and Bessel functions in the limit ν ≫ k_fr″ 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_fr″ = 200 that we encounter,

z_{ν}^{(3)} (k_{f} r^{″}) \approx 8.3 \times 10^{567}

and

z_{ν}^{(1)} (k_{f} r^{″}) \approx 3.05 \times 10^{- 574}

, while the maximum and minimum positive floating point numbers for double precision floating point format are 1.7976 × 10⁺³⁰⁸ and 2.225 × 10⁻³⁰⁸ 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. [20]) gives us a hint about the appropriate factors. It is observed that coefficients $D_{ν m}^{l M}$ and $D_{ν m}^{l N}$ which represent the coefficients of M and N waves in the scattered field in vacuum (Eq. 19 in Ref. [20]) are accompanied by factors of the form $z_{l}^{(1)} (k_{a} R_{1}) / z_{β}^{(1)} (k_{f} R_{2})$ and $z_{l}^{(1)} (k_{a} R_{1}) / ζ_{β}^{(1)} (k_{f} R_{2})$ respectively. Hence, what is important for the computation of the radiative energy transfer are not the terms $D_{ν m}^{l M}$ and $D_{ν m}^{l N}$ themselves but $D_{ν m}^{l M} z_{l}^{(1)} (k_{a} R_{1}) / z_{ν}^{(1)} (k_{f} R_{2})$ and $D_{ν m}^{l N} z_{l}^{(1)} (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. [20]) 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_{ν}^{(1)} (k_{f} R_{1}) A_{ν m, n m} / z_{n}^{(3)} (k_{f} R_{2})$ . Since A_νm,nm and B_νm,nm are related to β_νm,nm as given in Eqs. (17) and (18), it is sufficient if we build recursion relations for the normalized scalar translation coefficient $β_{ν m, n m}^{N}$ given by

β_{ν m, n m}^{N} = \frac{z_{ν}^{(1)} (k_{f} R_{1})}{z_{n}^{(3)} (k_{f} R_{2})} β_{ν m, n m} .

From Eq. (22) and Eq. (15), Eq. (16) it can be shown that the recursion relations to compute

β_{ν m, n m}^{N}

are given by :

\begin{array}{l} β_{ν m, n + 1 m}^{N} = \frac{1}{a_{n 0}^{+}} (- a_{n 0}^{-} \frac{z_{n + 1}^{(1)} (k_{f} R_{1})}{z_{n - 1}^{(1)} (k_{f} R_{1})} β_{ν m, n - 1 m}^{N} + a_{ν - 1, 0}^{+} \frac{z_{n + 1}^{(1)} (k_{f} R_{1})}{z_{n}^{(1)} (k_{f} R_{1})} \frac{z_{ν - 1}^{(3)} (k_{f} R_{2})}{z_{ν}^{(3)} (k_{f} R_{2})} β_{ν - 1 m, n m}^{N} + \\ a_{ν - 1, 0}^{-} \frac{z_{n + 1}^{(1)} (k_{f} R_{1})}{z_{n}^{(1)} (k_{f} R_{1})} \frac{z_{ν + 1}^{(3)} (k_{f} R_{2})}{z_{ν}^{(3)} (k_{f} R_{2})} β_{ν + 1 m, n m}^{N}), \end{array}

and

\begin{array}{l} β_{ν n + 1, n + 1 n + 1}^{N} = \frac{1}{b_{n n}^{+}} (b_{ν - 1, n}^{+} \frac{z_{n}^{(3)} (k_{f} R_{2})}{z_{n + 1}^{(3)} (k_{f} R_{2})} \frac{z_{ν}^{(1)} (k_{f} R_{1})}{z_{ν - 1}^{(1)} (k_{f} R_{1})} β_{ν - 1 n, n n}^{N} + \\ b_{ν + 1, n}^{-} \frac{z_{n}^{(3)} (k_{f} R_{2})}{z_{n + 1}^{(3)} (k_{f} R_{2})} \frac{z_{ν}^{(1)} (k_{f} R_{1})}{z_{ν + 1}^{(1)} (k_{f} R_{1})} β_{ν + 1 n, n n}) . \end{array}

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_fr″. The ratios have been computed by method of continued fractions [51]. The initialization for this recursion corresponding to Eq. (19) will be:

β_{ν 0, 00}^{N} = \frac{z_{0}^{(1)} (k_{f} R_{1})}{z_{ν}^{(3)} (k_{f} R_{2})} z_{ν}^{(3)} (k_{f} r^{″}) {(- 1)}^{ν} \sqrt{2 ν + 1} .

The spherical Hankel functions in Eq. (25), if computed individually, will lead to diverging numbers when ν ≫ k_fr″ as explained previously. To circumvent this,

β_{ν 0, 00}^{N}

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. [52].

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.

3. 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. [42] 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. (14). 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_fr″, the recurrence relations can be simplified further to make them more computationally viable.

Consider the form for the scalar translation coefficient β_νm,nm(k_fr″) given in Eq. (14). We note that the expression for β_νm,nm(k_fr″) depends on a summation over different orders p of the spherical Hankel function $z_{p}^{(3)} (k_{f} r^{″})$ . Since the magnitude of $z_{p}^{(3)} (k_{f} r^{″})$ increases exponentially for p ≫ k_fr″ (see Eq. (20)) it is possible to show that for p ≫ k_fr″ 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 :

β_{ν m, n m} (k_{f} r^{″}) = 4 π i^{2 ν} z_{n + ν}^{(3)} (k_{f} r^{″}) \sqrt{\frac{2 (n + ν) + 1}{4 π}} A (m, n, - m, ν, n + ν) .

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

\sqrt{2 p + 1}

vary relatively weakly with p, in comparison to

z_{p}^{(3)} (k_{f} r^{″})

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

ε \approx \frac{z_{n + ν - 2}^{(3)} (k_{f} r^{″})}{z_{n + ν}^{(3)} (k_{f} r^{″}) + z_{n + ν - 2}^{(3)} (k_{f} r^{″})} .

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

\frac{z_{p + 2}^{(3)} (k_{f} r^{″})}{z_{p}^{(3)} (k_{f} r^{″})} \approx {(\frac{2 p + 5}{k_{f} r^{″}})}^{2} .

The error from Eq. (27) then reduces to

ε \approx \frac{{(\frac{k_{f} r^{″}}{2 (n + ν) - 3})}^{2}}{{(\frac{k_{f} r^{″}}{2 (n + ν) - 3})}^{2} + 1} .

For 1 % error we get n + ν ≈ 5k_fr″. For our computations we have used n + ν ≥ 7k_fr″ as the criterion for employing the one-term approximation. From Fig. 2 it can be observed that when n + ν = 7k_fr″, $z_{n + ν}^{(3)} (k_{f} r^{″}) \approx 200 z_{n + ν - 2}^{(3)} (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. [43]) when the one-term approximation for β_νm,nm(k_fr″) is used beyond n + ν = pk_fr″ 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_fr″ as the criterion for employing the one-term approximation gives us an error of about 0.05 % in the spectral conductance.

Fig. 2 (a) Plot of ( $z_{n + ν}^{(3)} (k_{f} r^{″}) / z_{n + ν - 2}^{(3)} (k_{f} r^{″})$ ) as a function of n + ν for k_fr″ = 33 (arbitrarily chosen). The point n + ν = 7k_fr″, 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_fr″ are chosen as the criteria for employing the one-term approximation. From the plot, a criterion n + ν = 7k_fr″ 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₁ = 13.7 μm and R₂ = 40R₁ with minimum gap d/R₁ = 0.01.

Download Full Size | PDF

Since computation of Wigner 3j coefficient is computationally tedious, recursion relations can be derived for the approximated form for β_νm,nm. Employing $(\begin{matrix} n & ν & n + ν \\ - m & m & 0 \end{matrix}) = {(- 1)}^{n - ν} \sqrt{\frac{(2 n)! (2 ν)! (n + ν)! (n + ν)!}{(2 n + 2 ν + 1)! (n - m)! (n + m)! (ν - m)! (ν + m)!}}$ it is easy to show that the recursion relations turn out to be :

β_{ν m, n + 1 m} = β_{ν m, n m} \frac{n + 1}{n} \sqrt{\frac{(2 n + 3) (2 n + 1)}{{(n + 1)}^{2} - m^{2}}} \frac{n + ν + 1}{2 n + 2 ν + 1} \frac{z_{n + ν + 1}^{(3)} (k_{f} r^{″})}{z_{n + ν}^{(3)} (k_{f} r^{″})} .

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. (15). 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).

4. 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₁, R₂ (here, d = r″ − R₁ − R₂), the normalized translation coefficients $β_{n m, ν m}^{N}$ show a marked dependence on the radius ratio R₂/R₁. This behavior can be used to improve the efficiency of computation of near-field quantities considerably. In the one-term approximation, $β_{ν m, n m}^{N}$ varies with n and ν mainly through the ratio:

z_{n + ν}^{(3)} (k_{f} r^{″}) \frac{z_{ν}^{(1)} (k_{f} R_{1})}{z_{n}^{(3)} (k_{f} R_{2})} .

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₂/R₁ = 1, 3, 10, and 20, respectively, as a function of n and ν. The radius of the smaller sphere, R₁, and the minimum gap between the spheres, d, are taken to be 10 μm and 50 nm for all the cases. When R₂/R₁ = 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⁻⁶ 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₂/R₁ > 1 the distribution of the dominant terms in the matrix changes markedly. From Fig. 3(c) and 3(d), where R₂/R₁ ≫ 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

β_{ν m, n m}^{N}

) are compared in Table 1.

Fig. 3 Contour plots of the expression ${log}_{10} | z_{ν}^{(1)} (k_{f} R_{1}) z_{n + ν}^{(3)} (k_{f} r^{″}) / z_{n}^{(3)} (k_{f} R_{2}) |$ as a function of n and ν for two spheres with successive radius ratios R₂/R₁ = (a) 1, (b) 3, (c) 10, and (d) 20 with R₁ = 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

Download Full Size | PDF

Table 1. Table showing the error from the two approximations. Here, $G_{ω}^{Exact}$ 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⁻¹. 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₂ = R₁, and (2) R₂ = 10R₁. R₁ is kept a constant at 13.7 μm.

View Table

The variation of β_νm,nm with R₂/R₁ can be determined from the asymptotic behavior of the expression in Eq. (31), where k_fr″ = k_fR₁ + k_fR₂ + δ with δ = k_fd → 0 and n ≫ k_fR₁ and ν ≫ k_fR₂. From Ref. [48], the asymptotic form for $z_{ν}^{(1)} (k_{f} R_{1})$ for ν ≫ k_fR₁ is given by :

z_{ν}^{(1)} (k_{f} R_{1}) \approx \frac{exp [\sqrt{{(ν + 1 / 2)}^{2} - {(k_{f} R_{1})}^{2}} - (ν + 1 / 2) {cosh}^{- 1} [\frac{(ν + 1 / 2)}{k_{f} R_{1}}]]}{2 \sqrt{k_{f} R_{1} \sqrt{{(ν + 1 / 2)}^{2} - {(k_{f} R_{1})}^{2}}}} .

Similar asymptotic forms for

z_{n}^{(3)} (k_{f} R_{2})

and

z_{n + ν}^{(3)} (k_{f} r^{″})

can be obtained by replacing ν with n or n + ν. When n ≫ 1/2, ν ≫ 1/2, the magnitude of the expression in Eq. (31) is approximately:

\begin{matrix} exp [\sqrt{ν^{2} - {(k_{f} R_{1})}^{2}} + \sqrt{n^{2} - {(k_{f} R_{2})}^{2}} - \sqrt{{(n + ν)}^{2} - {(k_{f} R_{1} + k_{f} R_{2} + δ)}^{2}}] \times \\ exp [(n + ν) {cos}^{- 1} (\frac{n + ν}{k_{f} R_{1} + k_{f} R_{2} + δ}) - n {cosh}^{- 1} (\frac{n}{k_{f} R_{2}}) - ν {cosh}^{- 1} (\frac{ν}{k_{f} R_{1}})] . \end{matrix}

When δ reduces to zero, employing

{cosh}^{- 1} x = log (x + \sqrt{x^{2} - 1})

, taking n ≫ k_fR₁ and ν ≫ k_fR₂ and rearranging, the expression in Eq. (33) can be written as:

{(1 + \frac{n}{ν})}^{n} {(1 + \frac{ν}{n})}^{ν} {(1 + \frac{R_{2}}{R_{1}})}^{- ν} {(1 + \frac{R_{1}}{R_{2}})}^{- n} .

The expression in Eq. (34) can be shown to attain a maximum value of 1 when

\frac{n}{ν} = \frac{R_{2}}{R_{1}}

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 [- \frac{δ}{\sqrt{ν^{2} - {(k_{f} R_{1})}^{2}}} (\frac{ν^{2}}{k_{f} R_{1}} - k_{f} R_{1})] .

The above expression becomes negligibly small when

δ \sqrt{ν^{2} - {(k_{f} R_{1})}^{2}} ≫ k_{f} R_{1}

, i.e., if

ν ≫ \sqrt{{(\frac{R_{1}}{d})}^{2} + {(k_{f} R_{1})}^{2}},

From Eq. (35) the condition for n would be :

n ≫ \frac{R_{2}}{R_{1}} \sqrt{{(\frac{R_{1}}{d})}^{2} + {(k_{f} R_{1})}^{2}} .

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. [49] to be of the form:

N_{conv} = 2 \frac{R_{2}}{d} + e k_{f} \frac{r^{″}}{2},

This criterion was derived in Ref. [49] heuristically based on comparison with the convergence criterion for the near-field radiative transfer between planar bodies. Equations (37) and (38) 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).

5. 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. [20], modified suitably to account for the changes described in Secs. 2–4. The radiative conductance G (units WK⁻¹) between the two spheres shown in Fig. 1 is defined as:

G = lim_{T_{A} \to T_{B}} \frac{P (T_{A}, T_{B})}{| T_{A} - T_{B} |},

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 [43]. Here, we investigate the relevance of MPA for the case when R₂ ≫ R₁.

Fig. 4 (a) Application of MPA for the two spheres where the curved surfaces 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

Download Full Size | PDF

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

G^{MPA} (d, T) = \int_{0}^{R_{1}} h_{n f} (z) 2 π r d r + G_{c} (d, T),

where R₁ is the radius of the smaller sphere,

z = d + R_{1} - R_{2} - \sqrt{R_{1}^{2} - r^{2}} - \sqrt{R_{2}^{2} - r^{2}}

is the gap at distance r from the symmetry axis as shown in Fig. 4(a), h_nf (z) is the near-field contribution [43] 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 [54]:

G_{c} (d, T) = \frac{4 σ T^{3} (4 π R_{1}^{2})}{[(1 - ε) / ε] [1 + R_{1}^{2} / R_{2}^{2}] + 1 / F_{12} (d)},

where F₁₂ is the gap dependent view factor between the two spheres and can be approximated (to an accuracy of 1%) by [55]:

F_{12} (d) = \frac{1}{2} (1 - \sqrt{1 - \frac{1}{{(d / R_{2} + R_{1} / R_{2} + 1)}^{2}}}) .

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₁ = 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:

\frac{G_{c} (d_{1}, T)}{G_{c} (d_{2}, T)} = \frac{[(1 - ε) / ε] (1 + R_{1}^{2} / R_{2}^{2}) + 1 / F_{12} (d_{2})}{[(1 - ε) / ε] (1 + R_{1}^{2} / R_{2}^{2}) + 1 / F_{12} (d_{1})},

The comparison between G and G^MPA is shown in Fig. 5(a) for two spheres with R₁ = 13.7 μm and 2.5 μm and R₂ = 40R₁. The error between G and G^MPA is plotted in Fig. 5(b). For gaps d/R₁ < 0.1, MPA is observed to be able to predict the exact computed values of the conductance with errors less than 1%.

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₁ for two spheres with R₂/R₁ = 40. The study has been performed for R₁ = 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₁

Download Full Size | PDF

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 [43]), 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. (41) has to be appropriately modified to reflect spectral conductance. From proximity approximation theory, as d/R₁ → 0, the spectral conductance G_ω(d) for two spheres of unequal radii R₁ and R₂ is expected to scale as [56] :

G_{ω} (d) ~ \frac{R_{1} R_{2}}{R_{1} + R_{2}} \frac{1}{d} \approx \frac{R_{1}}{d} (for R_{2} ≫ R_{1}) .

This characteristic R₁/d behavior is observed for the resonant frequency contributions shown in Fig. 6(a) for d/R₁ ≲ 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).

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₁ for two spheres of radius R₁ = 2.5 μm, 13.7 μm and R₂ = 40R₁ (b) The % error between G_ω and $G_{ω}^{MPA}$ as a function of d/R₁

Download Full Size | PDF

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₁ for two spheres of radius R₁ = 2.5μm, 13.7μm and R₂ = 40R₁ (b) The % error between G_ω and $G_{ω}^{MPA}$ as a function of d/R₁

Download Full Size | PDF

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₁ ≳ 2 for R₁ = 13.7 μm, and d/R₁ ≳ 10 for R₁ = 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 with gap. As observed from Fig. 6(b) and Fig. 7(b) there is good agreement with the exact computed values of the conductance at such gaps. For intermediate gaps (2 ≲ d/R₁ ≲ 0.07, for R₁ = 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₁ ≲ 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₁, whereas at the non-resonant frequency the error is observed to be greater than 10% when R₁ = 2.5 μm. Despite such high errors at the non-resonant frequencies, MPA is successful in predicting the overall conductance when R₁ = 2.5 μm with error less than 1% for d/R₁ ≲ 0.1 as observed in Fig. 5. This apparent discrepancy has been analyzed in detail in Ref. [57] and has been attributed to decreased contribution from non-resonant frequencies as R₁ decreases.

We summarize the main contributions from this work:

The exponential behavior of the spherical Hankel function $z_{n}^{(3)} (k_{f} r^{″})$ when n ≥ k_fr″ 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_fr″ 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 and computationally less expensive compared to those for the exact forms of the translation coefficients.
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.
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₁ < 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.

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. [20]):

C_{n m}^{l M} + {\bar{u}}_{n} (R_{1}) \frac{z_{n}^{(1)} (k_{f} R_{1})}{z_{n}^{(3)} (k_{f} R_{1})} \sum_{ν = (m, 1)}^{N_{\max}} [D_{ν m}^{l M} A_{ν m, n m} (- k_{f} r^{″}) + D_{ν m}^{l N} B_{ν m, n m} (- k_{f} r^{″})] = p_{N}^{M} δ_{N l},

D_{n m}^{l M} + {\bar{u}}_{n} (R_{2}) \frac{z_{n}^{(1)} (k_{f} R_{2})}{z_{n}^{(3)} (k_{f} R_{2})} \sum_{ν = (m, 1)}^{N_{\max}} [C_{ν m}^{l M} A_{ν m, n m} (+ k_{f} r^{″}) + C_{ν m}^{l N} B_{ν m, n m} (+ k_{f} r^{″})] = 0,

C_{n m}^{l N} + {\bar{v}}_{n} (R_{1}) \frac{z_{n}^{(1)} (k_{f} R_{1})}{z_{n}^{(3)} (k_{f} R_{1})} \sum_{ν = (m, 1)}^{N_{\max}} [D_{ν m}^{l M} B_{ν m, n m} (- k_{f} r^{″}) + D_{ν m}^{l N} A_{ν m, n m} (- k_{f} r^{″})] = 0,

D_{n m}^{l N} + {\bar{v}}_{n} (R_{2}) \frac{z_{n}^{(1)} (k_{f} R_{2})}{z_{n}^{(3)} (k_{f} R_{2})} \sum_{ν = (m, 1)}^{N_{\max}} [C_{ν m}^{l M} B_{ν m, n m} (+ k_{f} r^{″}) + C_{ν m}^{l N} A_{ν m, n m} (+ k_{f} r^{″})] = 0,

where the symbol (m, 1) denotes the maximum of m and 1, ū_n(R₁) and v̄_n(R₁) are given by:

{\bar{u}}_{n} (R_{1}) = (\frac{k_{b} \frac{z_{n + 1}^{(1)} (k_{b} R_{1})}{z_{n}^{(1)} (k_{b} R_{1})} - k_{f} \frac{z_{n + 1}^{(1)} (k_{f} R_{1})}{z_{n}^{(1)} (k_{f} R_{1})}}{k_{b} \frac{z_{n + 1}^{(1)} (k_{b} R_{1})}{z_{n}^{(1)} (k_{b} R_{1})} - k_{f} \frac{z_{n + 1}^{(3)} (k_{f} R_{1})}{z_{n}^{(3)} (k_{f} R_{1})}}),

and

{\bar{v}}_{n} (R_{1}) = (\frac{n_{b} \frac{z_{n + 1}^{(1)} (k_{f} R_{1})}{z_{n}^{(1)} (k_{f} R_{1})} - \frac{z_{n + 1}^{(1)} (k_{b} R_{1})}{z_{n}^{(1)} (k_{b} R_{1})} + \frac{(n + 1) (1 - n_{b}^{2})}{k_{f} R_{1} n_{b}}}{n_{b} \frac{z_{n + 1}^{(1)} (k_{f} R_{1})}{z_{n}^{(1)} (k_{f} R_{1})} - \frac{z_{n + 1}^{(3)} (k_{b} R_{1})}{z_{n}^{(3)} (k_{b} R_{1})} + \frac{(n + 1) (1 - n_{b}^{2})}{k_{f} R_{1} n_{b}}}) .

ū_n(R₂) and v̄_n(R₂) have similar expressions. The normalizing factors for the coefficients

D_{ν m}^{l M}

and

D_{ν m}^{l N}

are

(z_{l}^{(1)} (k_{a} R_{1}) / z_{ν}^{(1)} (k_{f} R_{2}))

and

(z_{l}^{(1)} (k_{a} R_{1}) / ζ_{ν}^{(1)} (k_{f} R_{2}))

respectively. The equivalent normalizing factors for

C_{ν m}^{l M}

and

C_{ν m}^{l N}

are

(z_{l}^{(1)} (k_{b} R_{2}) / z_{n}^{(1)} (k_{f} R_{1}))

and

(z_{l}^{(1)} (k_{b} R_{2}) / ζ_{n}^{(1)} (k_{f} R_{1}))

respectively. Using these in Eqs. (46)–(49) and rearranging, the coupled linear equations reduce to:

\begin{array}{l} (\frac{z_{l}^{(1)} (k_{b} R_{2})}{z_{n}^{(1)} (k_{f} R_{1})} C_{n m}^{l M}) + {\bar{u}}_{n} (a) \frac{z_{l}^{(1)} (k_{b} R_{2})}{z_{l}^{(1)} (k_{a} R_{1})} \times \\ \sum_{ν = (m, 1)}^{N_{\max}} [\begin{matrix} (\frac{z_{l}^{(1)} (k_{a} R_{1})}{z_{ν}^{(1)} (k_{f} R_{2})} D_{ν m}^{l M}) (\frac{z_{ν}^{(1)} (k_{f} R_{2})}{z_{n}^{(3)} (k_{f} R_{1})} A_{ν m, n m} (- k_{f} D)) \\ + (\frac{z_{l}^{(1)} (k_{a} R_{1})}{ζ_{ν}^{(1)} (k_{f} R_{2})} D_{ν m}^{l N}) (\frac{ζ_{ν}^{(1)} (k_{f} R_{2})}{z_{n}^{(3)} (k_{f} R_{1})} B_{ν m, n m} (- k_{f} D)) \end{matrix}] = \frac{z_{l}^{(1)} (k_{b} R_{2})}{z_{n}^{(1)} (k_{f} R_{1})} p_{N}^{M} δ_{N l}, \end{array}

\begin{array}{l} (\frac{z_{l}^{(1)} (k_{a} R_{1})}{z_{n}^{(1)} (k_{f} R_{2})} D_{n m}^{l M}) + {\bar{u}}_{n} (R_{2}) \frac{z_{l}^{(1)} (k_{a} R_{1})}{z_{l}^{(1)} (k_{b} R_{2})} \times \\ \sum_{ν = (m, 1)}^{N_{\max}} [\begin{matrix} (\frac{z_{l}^{(1)} (k_{b} R_{2})}{z_{ν}^{(1)} (k_{f} R_{1})} C_{ν m}^{l M}) (\frac{z_{ν}^{(1)} (k_{f} R_{1})}{z_{n}^{(3)} (k_{f} R_{2})} A_{ν m, n m} (+ k_{f} D)) \\ + (\frac{z_{l}^{(1)} (k_{b} R_{2})}{ζ_{ν}^{(1)} (k_{f} R_{1})} C_{ν m}^{l N}) (\frac{ζ_{ν}^{(1)} (k_{f} R_{1})}{z_{n}^{(3)} (k_{f} R_{2})} B_{ν m, n m} (+ k_{f} D)) \end{matrix}] = 0, \end{array}

\begin{array}{l} (\frac{z_{l}^{(1)} (k_{b} R_{2})}{ζ_{n}^{(1)} (k_{f} R_{1})} C_{n m}^{l M}) + {\bar{v}}_{n} (R_{1}) \frac{z_{n}^{(1)} (k_{f} R_{1})}{ζ_{n}^{(1)} (k_{f} R_{1})} \frac{z_{l}^{(1)} (k_{b} R_{2})}{z_{l}^{(1)} (k_{a} R_{1})} \times \\ \sum_{ν = (m, 1)}^{N_{\max}} [\begin{matrix} (\frac{z_{l}^{(1)} (k_{a} R_{1})}{z_{ν}^{(1)} (k_{f} R_{2})} D_{ν m}^{l M}) (\frac{z_{ν}^{(1)} (k_{f} R_{2})}{z_{n}^{(3)} (k_{f} R_{1})} B_{ν m, n m} (- k_{f} D)) \\ + (\frac{z_{l}^{(1)} (k_{a} R_{1})}{ζ_{ν}^{(1)} (k_{f} R_{2})} D_{ν m}^{l N}) (\frac{ζ_{ν}^{(1)} (k_{f} R_{2})}{z_{n}^{(3)} (k_{f} R_{1})} A_{ν m, n m} (- k_{f} D)) \end{matrix}] = 0, \end{array}

\begin{array}{l} (\frac{z_{l}^{(1)} (k_{a} R_{1})}{ζ_{n}^{(1)} (k_{f} R_{2})} D_{n m}^{l M}) + {\bar{v}}_{n} (R_{2}) \frac{z_{n}^{(1)} (k_{f} R_{2})}{ζ_{n}^{(1)} (k_{f} R_{2})} \frac{z_{l}^{(1)} (k_{a} R_{1})}{z_{l}^{(1)} (k_{b} R_{2})} \times \\ \sum_{ν = (m, 1)}^{N_{\max}} [\begin{matrix} (\frac{z_{l}^{(1)} (k_{b} R_{2})}{z_{ν}^{(1)} (k_{f} R_{1})} C_{ν m}^{l M}) (\frac{z_{ν}^{(1)} (k_{f} R_{1})}{z_{n}^{(3)} (k_{f} R_{2})} B_{ν m, n m} (+ k_{f} D)) \\ + (\frac{z_{l}^{(1)} (k_{b} R_{2})}{ζ_{ν}^{(1)} (k_{f} R_{1})} C_{ν m}^{l N}) (\frac{ζ_{ν}^{(1)} (k_{f} R_{1})}{z_{n}^{(3)} (k_{f} R_{2})} A_{ν m, n m} (+ k_{f} D)) \end{matrix}] = 0 . \end{array}

Thus there are two possible normalizing factors for A_νm,nm(+k_fD):

(ζ_{ν}^{(1)} (k_{f} R_{1}) / z_{n}^{(3)} (k_{f} R_{2}))

and

(z_{ν}^{(1)} (k_{f} R_{1}) / z_{n}^{(3)} (k_{f} R_{2}))

. Since the function

ζ_{ν}^{(1)} (k_{f} R_{1})

behaves similar to

z_{ν}^{(1)} (k_{f} R_{1})

for ν ≫ k_fR₁ only one of them (the latter) has been chosen as a representative form for discussion in Eq. (22).

Acknowledgments

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

References and links

1. A. Narayanaswamy and G. Chen, “Surface modes for near field thermophotovoltaics,” Appl. Phys. Lett. 82, 3544–3546 (2003). [CrossRef]

2. S. Basu, Z. M. Zhang, and C. J. Fu, “Review of near-field thermal radiation and its application to energy conversion,” Int. J. Energy Res. 33(13), 1203–1232 (2009). [CrossRef]

3. O. Ilic, M. Jablan, J. D. Joannopoulos, I. Celanovic, and M. Soljačić, “Overcoming the black body limit in plasmonic and graphene near-field thermophotovoltaic systems,” Opt. Express 20(103), A366–A384 (2012). [CrossRef] [PubMed]

4. M. Laroche, R. Carminati, and J.-J. Greffet, “Near-field thermophotovoltaic energy conversion,” J. Appl. Phys. 100(6), 063704 (2006). [CrossRef]

5. C. Otey, W.-T. Lau, and S. Fan, “Thermal rectification through vacuum,” Phys. Rev. Lett. 104(15), 154301 (2010). [CrossRef] [PubMed]

6. B. J. Lee, Y.-B. Chen, and Z. M. Zhang, “Confinement of infrared radiation to nanometer scales through metallic slit arrays,” J. Quant. Spectrosc. Radiat. Transfer 109(4), 608–619 (2008). [CrossRef]

7. W. A. Challener, C. Peng, A. V. Itagi, D. Karns, W. Peng, Y. Peng, X. Yang, X. Zhu, N. J. Gokemeijer, and Y. T. Hsia, “Heat-assisted magnetic recording by a near-field transducer with efficient optical energy transfer,” Nat. Photon. 3(4), 220–224 (2009). [CrossRef]

8. Y. D. Wilde, F. Formanek, R. Carminati, B. Gralak, P. A. Lemoine, K. Joulain, J. P. Mulet, Y. Chen, and J.-J. Greffet, “Thermal radiation scanning tunnelling microscopy,” Nature 444(7120), 740–743 (2006). [CrossRef] [PubMed]

9. B. Guha, C. Otey, C. B. Poitras, S. Fan, and M. Lipson, “Near-field radiative cooling of nanostructuresm,” Nano Lett. 12(9), 4546–4550 (2012). [CrossRef] [PubMed]

10. S. M. Rytov, “Theory of electric fluctuations and thermal radiation,” Air Force Cambridge Research Center, Bedford, MA (1959).

11. D. Polder and M. V. Hove, “Theory of radiative heat transfer between closely spaced bodies,” Phys. Rev. B 4, 3303–3314 (1971). [CrossRef]

12. J. J. Loomis and H. J. Maris, “Theory of heat transfer by evanescent electromagnetic waves,” Phys. Rev. B 50, 18517–18524 (1994). [CrossRef]

13. E. G. Cravalho, C. L. Tien, and R. P. Caren, “Effect of small spacings on radiative transfer between two dielectrics,” J. Heat Transfer 89, 351–358 (1967). [CrossRef]

14. J.-P. Mulet, K. Joulain, R. Carminati, and J.-J. Greffet, “Enhanced radiative transfer at nanometric distances,” Microscale Thermophys. Eng. 6, 209–222 (2002). [CrossRef]

15. J. B. Pendry, “Radiative exchange of heat between nanostructures,” J. Phys.: Condens. Matter 11, 6621 (1999).

16. C. H. Park, H. A. Haus, and M. S. Weinberg, “Proximity-enhanced thermal radiation,” J. Phys. D: Appl. Phys 35(21), 2857–2863 (2002). [CrossRef]

17. A. I. Volokitin and B. N. J. Persson, “Radiative heat transfer between nanostructures,” Phys. Rev. B 63, 205404 (2001). [CrossRef]

18. G. Domingues, S. Volz, K. Joulain, and J.-J. Greffet, “Heat transfer between two nanoparticles through near-field interaction,” Phys. Rev. Lett. 94(8), 085901 (2005). [CrossRef] [PubMed]

19. J.-P. Mulet, K. Joulain, R. Carminati, and J.-J. Greffet, “Nanoscale radiative heat transfer between a small particle and a plane surface,” Appl. Phys. Lett. 78, 2931–2933 (2001). [CrossRef]

20. A. Narayanaswamy and G. Chen, “Thermal near-field radiative transfer between two spheres,” Phys. Rev. B 77(7), 075125 (2008). [CrossRef]

21. C. Otey and S. Fan, “Numerically exact calculation of electromagnetic heat transfer between a dielectric sphere and plate,” Phys. Rev. B 84, 245431 (2011). [CrossRef]

22. M. Krüger, T. Emig, and M. Kardar, “Nonequilibrium electromagnetic fluctuations: Heat transfer and interactions,” Phys. Rev. Lett. 106(21), 210404 (2011). [CrossRef] [PubMed]

23. A. W. Rodriguez, M. T. H. Reid, and S. G. Johnson, “Fluctuating-surface-current formulation of radiative heat transfer for arbitrary geometries,” Phys. Rev. B 86(22), 220302 (2012). [CrossRef]

24. R. Messina and M. Antezza, “Scattering-matrix approach to Casimir-Lifshitz force and heat transfer out of thermal equilibrium between arbitrary bodies,” Phys. Rev. A 84(4), 042102 (2011). [CrossRef]

25. A. Narayanaswamy and Y. Zheng, “A Green’s function formalism of energy and momentum transfer in fluctuational electrodynamics,” J. Quant. Spectrosc. Radiat. Transfer 132, 12–21 (2014). [CrossRef]

26. P. Ben-Abdallah, S. Biehs, and K. Joulain, “Many-body radiative heat transfer theory,” Phys. Rev. Lett. , 107(11), 114301 (2011). [CrossRef] [PubMed]

27. S. Shen, A. Narayanaswamy, and G. Chen, “Surface phonon polaritons mediated energy transfer between nanoscale gaps,” Nano Lett. 9(8), 2909–2913 (2009). [CrossRef] [PubMed]

28. E. Rousseau, A. Siria, G. Jourdan, S. Volz, F. Comin, J. Chevrier, and J.-J. Greffet, “Radiative heat transfer at the nanoscale,” Nat. Photon. 3(9), 514–517 (2009). [CrossRef]

29. P. J. van Zwol, L. Ranno, and J. Chevrier, “Tuning near field radiative heat flux through surface excitations with a metal insulator transition,” Phys. Rev. Lett. 108(23), 234301 (2012). [CrossRef] [PubMed]

30. S. Shen, A. Mavrokefalos, P. Sambegoro, and G. Chen, “Nanoscale thermal radiation between two gold surfaces,” Appl. Phys. Lett. 100(23), 233114 (2012). [CrossRef]

31. A. Narayanaswamy, S. Shen, and G. Chen, “Near-field radiative heat transfer between a sphere and a substrate,” Phys. Rev. B 78(11), 115303 (2008). [CrossRef]

32. W. M. Hirsch, A. Kraft, M. T. Hirsch, J. Parisi, and A. Kittel, “Heat transfer in ultrahigh vacuum scanning thermal microscopy,” J. Vac. Sci. Technol. A 17(4), 1205–1210 (1999). [CrossRef]

33. A. Kittel, W. Müller-Hirsch, J. Parisi, S. A. Biehs, D. Reddig, and M. Holthaus, “Near-field heat transfer in a scanning thermal microscope,” Phys. Rev. Lett. 95(22), 224301 (2005). [CrossRef] [PubMed]

34. L. Worbes, D. Hellmann, and A. Kittel, “Enhanced near-field heat flow of a monolayer dielectric island,” Phys. Rev. Lett. 110(13), 134302 (2013). [CrossRef] [PubMed]

35. A. D. Yaghjian, “Electric dyadic Green’s functions in the source region,” Proc. IEEE 68(2), 248–263 (1980). [CrossRef]

36. R. E. Collin, Field Theory of Guided Waves (IEEE Press, 1991), Vol. 2.

37. L. Landau, E. M. Lifšic, J. B. Sykes, J. S. Bell, M. J. Kearsley, and L. Pitaevskii, Electrodynamics of Continuous Media (Pergamon, 1960).

38. R. Kubo, “The fluctuation-dissipation theorem,” Rep. Prog. Phys. 29(1), 255 (1966). [CrossRef]

39. U. M. B. Marconi, A. Puglisi, L. Rondoni, and A. Vulpiani, “Fluctuation–dissipation: response theory in statistical physics,” Phys. Rep. 461(4), 111–195 (2008). [CrossRef]

40. J. H. Bruning and Y. T. Lo, “Multiple scattering by spheres,” Antenna Laboratory Report No.69(5) (1969).

41. W. C. Chew, “Derivation of the vector addition theorem,” Microwave Opt. Technol. Lett. 3, 256–260 (1990). [CrossRef]

42. W. C. Chew, “Efficient ways to compute the vector addition theorem,” J. Electromagn. Wave 7, 651–665 (1993). [CrossRef]

43. K. Sasihithlu and A. Narayanaswamy, “Proximity effects in radiative heat transfer,” Phys. Rev. B 83(16), 161406 (2011). [CrossRef]

44. J. Błocki, J. Randrup, W. J. Światecki, and C. F. Tsang, “Proximity forces,” Ann. Phys. 105(2), 427–462 (1977). [CrossRef]

45. S. K. Lamoreaux, “Demonstration of the Casimir force in the 0.6 to 6μm range,” Phys. Rev. Lett. 78, 5–8 (1997). [CrossRef]

46. H. Gies and K. Klingmüller, “Casimir effect for curved geometries: Proximity-force-approximation validity limits,” Phys. Rev. Lett. 96(22), 220401 (2006). [CrossRef] [PubMed]

47. V. Golyk, M. Krüger, A. P. McCauley, and M. Kardar, “Small distance expansion for radiative heat transfer between curved objects,” Europhys. Lett. 101(3), 34002 (2013). [CrossRef]

48. M. Abramowitz and I. Stegun, Handbook of Mathematical Functions: With Formulas, Graphs, and Mathematical Tables (Dover, 1965).

49. K. Sasihithlu and A. Narayanaswamy, “Convergence of vector spherical wave expansion method applied to near-field radiative transfer,” Opt. Express 19(S4), A772–A785 (2011). [CrossRef] [PubMed]

50. W. C. Chew, “Recurrence relations for three-dimensional scalar addition theorem,” J. Electromagn. Wave 6, 133–142 (1992). [CrossRef]

51. A. Cuyt, Handbook of Continued Fractions for Special Functions (Springer, 2008).

52. E. J. Rothwell, “Computation of the logarithm of bessel functions of complex argument,” Commun. Numer. Methods Eng. 21(10), 597–605 (2005). [CrossRef]

53. W. C. Chew, “A derivation of the vector addition theorem,” Microwave Opt. Tech. Lett. 3(7), 256–260 (1990). [CrossRef]

54. R. Siegel and J. R. Howell, Thermal Radiation Heat Transfer (Taylor and Francis, 2002).

55. N. H. Juul, “Investigation of approximate methods for calculation of the diffuse radiation configuration view factor between two spheres,” Lett. Heat Mass Transfer 3(6), 513–521 (1976). [CrossRef]

56. B. E. Sernelius, Surface Modes in Physics (Wiley-Vch, 2011).

57. N. Gu, K. Sasihithlu, Y. Zheng, and A. Narayansawamy, “Proximity approximation and radiative transfer between sub-wavelength spheres and planar surfaces,” under review (2013).

Case	$G_{ω}^{Exact}$	$G_{ω}^{A_{1}}$	Error A1	$G_{ω}^{A_{2}}$	Error A2
R₂ = R₁	132.1995	133.799	1.21	132.1745	1.9 × 10⁻²
R₂ = 10R₁	242.6136	241.0964	0.625	242.6132	1.6 × 10⁻⁴

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

Abstract

1. Introduction

2. Translation addition theorems

3. One term approximation

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

5. Validity of modified proximity approximation

Appendix: Normalizing factors for the vector translation coefficients A_νm,nm and B_νm,nm

Acknowledgments

References and links

Cited By

Figures (7)

Tables (1)

Equations (57)

Optics Express

Abstract

1. Introduction

2. Translation addition theorems

3. One term approximation

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

5. Validity of modified proximity approximation

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

Acknowledgments

References and links

Cited By

Figures (7)

Tables (1)

Equations (57)

Optics Express

Appendix: Normalizing factors for the vector translation coefficients A_νm,nm and B_νm,nm