Discrete dipole approximation for the study of radiation dynamics in a magnetodielectric environment

We develop a general computational approach, based on the discrete dipole approximation, for the study of radiation dynamics near or inside an object with arbitrary linear dielectric permittivity, and magnetic permeability tensors. Our method can account for dispersion and losses and provides insight on the role of local-field corrections in discrete magnetodielectric structures. We illustrate our method in the case of a source inside a magneto-dielectric, isotropic sphere for which the spontaneous emission rate of a source can be computed analytically. We show that our approach is in excellent agreement with the exact result, providing an approach capable of handling both the electric and magnetic response of advanced metamaterials. © 2010 Optical Society of America OCIS codes: (160.3918) Metamaterials; (290.0290) Scattering. References and links 1. E. Purcell, “Spontaneous emission probabilities at radio frequencies,” Phys. Rev. 69, 681 (1946). 2. K. Drexhage, “Influence of a dielectric interface on fluorescence decay time,” J. Lumin. 1, 693 (1970). 3. S. Haroche, “Cavity quantum electrodynamics,” in Fundamental Systems in Quantum Optics, J. Dalibard, J. M. Raimond, and J. Zinn-Justin, eds. (North-Holland, Amsterdam, 1992), pp. 767–940. 4. N. Litchinitser and V. Shalaev, “Photonic metamaterials,” Laser Phys. Lett. 5, 411–420 (2008). 5. E. M. Purcell and C. R. Pennypacker, “Scattering and absorption of light by nonspherical dielectric grains,” Astrophys. J. 186, 705–714 (1973). 6. B. T. Draine, “The discrete-dipole approximation and its application to interstellar graphite grains,” Astrophys. J. 333, 848–872 (1988). 7. B. T. Draine and P. J. Flatau, “Discrete-dipole approximation for scattering calculations,” J. Opt. Soc. Am. A 11, 1491–1499 (1994). 8. A. Rahmani, P. Chaumet, and G. Bryant, “Coupled dipole method with an exact long-wavelength limit and improved accuracy at finite frequencies,” Opt. Lett. 27, 2118–2120 (2002). 9. A. Rahmani, P. C. Chaumet, and G. W. Bryant, “On the Importance of Local-Field Corrections for Polarizable Particles on a Finite Lattice: Application to the Discrete Dipole Approximation,” Astrophys. J. 607, 873–878 (2004). 10. M. A. Yurkin, V. P. Maltsev, and A. G. Hoekstra, “The discrete dipole approximation for simulation of light scattering by particles much larger than the wavelength,” J. Quant. Spect. Rad. Transf. 106, 546–557 (2007). 11. M. A. Yurkin and A. G. Hoekstra, “The discrete dipole approximation: An overview and recent developments,” J. Quant. Spect. Rad. Transf. 106, 558–589 (2007). 12. P. C. Chaumet, A. Rahmani, and G. W. Bryant, “Generalization of the coupled dipole method to periodic structure,” Phys. Rev. B 67, 165404 (2003). #125208 $15.00 USD Received 9 Mar 2010; accepted 30 Mar 2010; published 7 Apr 2010 (C) 2010 OSA 12 April 2010 / Vol. 18, No. 8 / OPTICS EXPRESS 8499 13. B. T. Draine and P. J. Flatau, “The discrete dipole approximation for periodic targets: Theory and tests,” J. Opt. Soc. Am. A 25 (2008). 14. B. T. Draine and J. C. Weingartner, “Radiative Torques on Interstellar Grains: I. Superthermal Spinup,” Astrophys. J. 470, 551–565 (1996). 15. A.G. Hoekstra, M. Frijlink, L.B.F.M. Waters, and P.M.A. Sloot, “Radiation forces in the discrete-dipole approximation,” J. Opt. Soc. Am. A 18, 1944–1953 (2001). 16. M. Nieto-Vesperinas, P. C. Chaumet, and A. Rahmani, “Near-field photonic forces,” Phil. Trans. Roy. Soc. Lond. A 362, 719–737 (2004). 17. P. C. Chaumet, A. Rahmani, and M. Nieto-Vesperinas, “Optical trapping and manipulation of nano-object with an apertureless probe,” Phys. Rev. Lett. 88, 123601 (2002). 18. A. Rahmani and P. C. Chaumet, “Optical Trapping near a Photonic Crystal,” Opt. Express 14, 6353–6358 (2006). http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-14-13-6353 19. A. Rahmani, P. C. Chaumet, and F. de Fornel, “Environment-induced modification of spontaneous emission: Single-molecule near-field probe,” Phys. Rev. A 63, 023819 (2001). 20. A. Rahmani and G. W. Bryant, “Spontaneous emission in microcavity electrodynamics,” Phys. Rev. A 65, 033817 (2002). 21. F. Bordas, N. Louvion, S. Callard, P. C. Chaumet, and A. Rahmani, “Coupled dipole method for radiation dynamics in finite photonic crystal structures,” Phys. Rev. E 73, 056601 (2006). 22. A. Lakhtakia, “Strong and weak forms of the method of moments and the coupled dipole method for scattering of time-harmonic electromagnetic fields,” Int. J. Mod. Phys. C 3, 583–603 (1992). 23. Y. You, G. W. Kattawar, P.-W. Zhai, and p. Yang, “Zero-backscatter cloak for aspherical particles using a generalized dda formalism,” Opt. Express 16, 2068–2079 (2008), http://www.opticsinfobase.org/oe/ abstract.cfm?URI=oe-16-3-2068cfm?URI=oe-16-3-2068 24. P. C. Chaumet and A. Rahmani, “Coupled-dipole method for magnetic and negative refraction materials,” J. Quant. Spect. Rad. Transf. 110, 22–29 (2009). 25. P. C. Chaumet and A. Rahmani, “Electromagnetic force and torque on magnetic and negative-index scatterers,” Opt. Express 17, 2224–2234 (2009), http://www.opticsinfobase.org/oe/abstract.cfm?URI=


Introduction
In 1946 Purcell pointed out that the dynamics of an electromagnetic source is not an intrinsic feature of the source but is, rather, a consequence of the coupling between the source and its environment [1].This idea, however, only became widely accepted in the late 1960s, in large part owing to Drexhage's pioneering work on emitters near mirrors [2].We now know that the emission rate of an emitter in the presence of matter will, in general, be different from the emission rate for the same emitter in vacuum [3].Traditionally, non-magnetic materials are used to control the radiation dynamics of electromagnetic sources, however, recent advances in metamaterial research have made it apparent that the ability to tailor both the electric and magnetic response of a material leads to a much larger wealth of electromagnetic behavior [4].
The first important step in the study of radiation dynamics in active metamaterials is to develop a general computational approach for radiation dynamics in magneto-dielectric (MD) media, i.e., media with both electric and magnetic properties.Moreover, achieving any significant control over the electric and magnetic properties of an electromagnetic structure would most likely require a composite system, with a discrete structure.Therefore it is also essential to understand how the discreteness of the structure may affect the radiation dynamics of a source.In this article we present an approach that fulfills these objectives.
Our approach is based on the discrete dipole approximation (DDA) (also called the coupled dipole method (CDM)) [5][6][7][8][9][10] (see [11] for review).The method was originally developed to study the scattering of light by dielectric interstellar dust particles.However, over the years the DDA has been used to study absorbing dielectrics and metals, and the method has been generalized to deal with periodic structures [12,13], optical forces [14][15][16], optical trapping [17,18], and radiation dynamics [19][20][21].While the DDA has been used traditionally to study the scattering of light by non-magnetic structures, the method has also been formulated to address electromagnetic scattering [22][23][24], and optical forces [25], for magnetic and negativeindex materials.
In this article, we present a further generalization of the DDA that can be used to compute the radiation dynamics of a source in an arbitrary environment with linear dielectric and magnetic responses.The generalization provides an approach capable of describing complex, advanced metamaterial structures.The example of a source in a spherical magnetodielectric cavity is presented to illustrate the accuracy of the generalized approach.Furthermore, this example highlights the importance of including magnetic effects.In particular, we will show that, in the continuous case, strong enhancement of the response of an electric source arises when the cavity has a magnetic response, and vice-versa, whereas, in the discrete case, the effect is mitigated by local-field enhancements.

DDA for Radiation Dynamics in an Arbitrary Magnetodielectric cavity
Consider an electromagnetic source located at position r 0 in the presence of a cavity (we shall use the term cavity as a generic term for any object near or containing the source).Let k = ω/c be the wave number where ω is the angular frequency and c is the speed of light in vacuum.Let us assume that the source undergoes an electric dipole transition with unit dipole moment p.As a shorthand we will refer to this situation as an electric source.In the weak coupling regime, the rate at which the source loses energy, radiatively and non-radiatively, in the presence of the cavity (Γ cav ), normalized to the free-space case (Γ ∞ ), is given by [20]: where Ḡcav is the sum of the electric and magnetic field-susceptibility tensors (FSTs) associated with the cavity.Therefore, the radiation dynamics of the source can be studied provided one knows the FST associated with the environment of the source.To this end, we need to establish a procedure to derive the FST, in the general case, using the DDA.
Let us consider an object with linear dielectric permittivity tensor ε(r, ω) and linear magnetic permeability tensor μ(r, ω).We discretize the object into N subunits over a cubic lattice with spacing d.Each subunit is assumed to be small enough to satisfy the dipole approximation.Omitting the dependence on ω henceforth, the optical response of subunit i is described by electric and magnetic polarizability tensors: with where Ī is the identity tensor.Let Ḡ ee 0 (r, r ) and Ḡ me 0 (r, r ) be free-space FSTs that give the electric and magnetic fields at r due to an electric source at r [26].Similarly, we introduce Ḡ em 0 (r, r ) and Ḡ mm 0 (r, r ) in the case of a magnetic source.The FSTs representing the linear response of the electromagnetic fields at subunit i due to an electric source at 0 , in the presence of the cavity, can be written in a self-consistent form: From the expressions above we can find Ḡ ee (r i , r 0 ) and Ḡ me (r i , r 0 ) by solving a linear system using, for instance, the techniques described in [27].The FST giving the electric field scattered by the object back to the source can now be expressed as: Inserting Ḡ ee cav (r 0 , r 0 ) in Eq. ( 1) yields the decay rate for an electric source in the presence of the object.A similar treatment would lead to Ḡ mm cav (r 0 , r 0 ) and to the decay rate for a magnetic source.

Local-field factor
We can illustrate our approach in the case of a source at the center of a lossless, homogeneous spherical magnetodielectric cavity with radius a.For this configuration, Ḡcav can be expressed analytically and the decay rate can be calculated following Chew's approach [28], allowing us to compare the DDA against the exact result.However, we should first note one essential difference between the two approaches.The analytical approach assumes that the cavity is made of a continuous medium.However, in the DDA treatment, and most likely in any actual metamaterial structure, the cavity is a composite medium with a discrete structure.This implies that we need to consider the influence of local-field corrections on the radiation dynamics of the source.Indeed, in a discrete structure, the radiation dynamics of a source may be altered by two types of boundary effects.The first one, which we have focused on so far, is the usual Purcell effect which relates to the shape of the outer boundary of the cavity.The second boundary effect corresponds to the fact that the source is placed between the elements forming the composite material (in our case the DDA subunits), rather than in a continuous medium.As a result, the electromagnetic field of the source will interact with two sets of boundaries, the "local" ones between the vacuum region where the source is and the subunits in the immediate neighborhood of the source, and the "global" ones which define the shape of the cavity.In the DDA, because we discretize the cavity, both contributions are actually included in Eq. (7).Accordingly, to compare the DDA result to an exact derivation of the decay rate in the case of a spherical cavity made of a continuous medium, we first need to compensate for the local-field effect (discreteness) which means isolating the contribution to the decay rate due to the discrete lattice.
As we showed in Ref. [29] in the case of an isotropic, lossless non-magnetic medium, if the source is introduced in the cavity without disturbing the DDA lattice, the local-field factor can be derived analytically as fast-converging lattice sums.In the present case, for a magnetodielectric medium, because the local field factor is derived in the quasistatic approximation, electric and magnetic effects are decoupled and, as a result, for an electric (magnetic) source the local-field factor only depends on the permittivity (permeability) of the medium.In other words, this means that for an electric source the expressions derived in Ref. [29] still hold, and for a magnetic source, all one has to do is replace ε by μ in the same expressions.
Of course, here we are assuming that the discrete forming the object are small enough that the induced polarization on each element is uniform.If the discrete elements had more complex geometry and scattering properties (e.g., the element are too large to satisfy the dipole approximation), the local-field correction would have to be computed by taking into account the actual electromagnetic properties of the discrete elements.However, such issues are beyond the scope of this work.Once the local-field factor is known, the DDA result can be normalized by this factor to produce the decay rate associated with a continous medium.The local-field factor depends on the position of the source within the DDA lattice.For instance, if we place the electric source at the center of a homogeneous spherical cavity in such a way that it lies at a reduced position ( 1 / 2 , 1 / 2 , 1 / 2 ) within the lattice the local-field factor will be (ε + 2)/3 [29].Note that this local-field factor does not depend on the lattice spacing provided that the spacing remains smaller than the wavelength inside the material and the material is lossless.If the material is lossy, the DDA can still be used to compute decay rates but the localfield effect will have a non-radiative component that will depend on the lattice spacing [20].

Continuous picture
To illustrate how the material properties of the cavity affect the spontaneous emission rate of the source we consider an electric source at the center of a lossless, homogeneous, isotropic spherical cavity with radius a and a fixed refractive index n = √ εμ = 2.However, while the refractive index is fixed, we consider different "apportionments" of the refractive index between the electric and magnetic responses of the cavity.We calculate, using both our DDA approach with cubic lattice spacing a/30 and Chew's theory [28], the normalized decay rate for the source as a function of the size of the sphere expressed in terms of the parameter ka for three values of (ε, μ): (4,1), (2,2), and (1,4).In Fig. 1(a), the analytical results are plotted as solid lines while the DDA results are plotted as open circles.The agreement of the DDA calculation with the exact result is excellent.For the three cases we observe the usual oscillations of the decay rate with the cavity size (interferences between an oscillating dipole and the reflected field).However, the most interesting feature emerges when we examine the magnitude of the decay rate, as we decrease (increase) the permittivity (permeability).Keeping in mind that the curves in Fig. 1(a) pertain to cavities with the same refractive index, we see that, over the range of cavity sizes considered here, the largest emission rate enhancement for an electric source is more than 5 times larger for ε = 1 and μ = 4 than for ε = 4 and μ = 1.The difference in the enhancement is a consequence of the fact that ε and μ do not appear on the same footing in the expression of the decay rate.This can be best illustrated in the simpler case of a source in an infinite medium with refractive index n.In the non-magnetic case (μ = 1) the decay rates for an electric source in the medium and in vacuum are related by Γ med = nΓ vac = √ εΓ vac .In the general case, however, we have Γ med = μnΓ vac = μ √ εμΓ vac [30].The dependence of Γ med on μ means that for a fixed value of the refractive index, the decay rate for an electric source in an infinite medium scales linearly with the magnetic permeability.On the other hand, the decay rate would scale as the inverse of the dielectric permittivity.The case of a magnetic source (i.e., a magnetic dipole transition) leads to a similar result with the interchange of ε and μ.In the case of a spherical cavity, the dependence of the decay rate on the permittivity (or the permeability) will depend on the size of the cavity, however, for a given refractive index, the largest enhancement of the emission rate for an electric source in a continuous spherical cavity is always achieved for μ = n 2 due to the dominant role played by the permeability in the expression of the decay rate.Thus if one is able to tailor the electric and magnetic response of the cavity, for a given refractive index of the cavity, it appears to be more beneficial to enhance the magnetic properties of the cavity in order to enhance the radiation dynamics of an electric source and vice-versa.

Discrete picture
If we do not normalize the DDA result with the local-field factor we get the decay rate for the discrete spherical cavity.The results are plotted as solid lines in Fig. 1(b).Note that the CDM result for (ε, μ)=(1,4) is the same as for the continuous cavity because in that case the local-field factor is 1 since for an electric source the local-field factor depends on ε only.For a discrete structure, the difference in rate enhancements for the three sets of material constants is less pronounced.Although a large μ yields an enhancement of the decay rate for an electric source, a large ε contributes to the enhancement through the local-field factor squared.This example highlights the interplay between Purcell and local-field effects in a discrete MD cavity.It also suggests that a study of radiation dynamics in electromagnetic metamaterials would most likely need to consider the discrete nature of the structure.

Conclusion
In conclusion we have developed a general formalism to study radiation dynamics in complex structures with linear permittivity and permeability tensors.We have validated our approach against the analytical result for a homogeneous magneto-dielectric spherical cavity.Furthermore, we have highlighted the role of local-field corrections in discrete structures.Although in our case each discretization cell supports an electric and a magnetic dipole resonance, in an actual metamaterial structure the form of the local field correction will depend on the specific geometry and material composition of each discrete element forming the metamaterial.While we focused here on the dissipative aspect (decay rate) of the interaction between the source and its environment, the dispersive aspect (frequency shift or anomalous Lamb shift) can also be derived once the field-susceptibility is known [20].Finally, we have seen that even a simple configuration such as an isotropic, homogeneous MD sphere can affect the dynamics of a source in an interesting way.This highlights the fact that it is essential to include the magnetic response of the material (the refractive index is not enough) in the description of the coupling between a source and a MD cavity.

Fig. 1 .
Fig. 1.Normalized emission rate as a function of the size parameter ka for an electric dipole at the center of a spherical cavity with radius a.(a) Solid lines: analytical calculation for a continuous, homogeneous cavity; open symbols: DDA calculation for a continuous, homogeneous cavity.The corresponding material constants (ε, μ) are given in the inset.(b)DDA calculation for a discrete, homogeneous cavity.The corresponding material constants (ε, μ) are given in the inset.The source is located at the reduced position ( 1 / 2 , 1 / 2 , 1 / 2 ) within a discretization cell.