Topology optimized all-dielectric cloak: design, performances and modal picture of the invisibility effect

Benjamin Vial; Yang Hao

doi:10.1364/OE.23.023551

1. Introduction

In recent years, the advent of Transformation Optics [1, 2] (TO) has paved the way for controlling light in an unprecedented way. The most striking example that has caught popular attention is the so-called invisibility cloak, a device that can bend light around an object to make it seem invisible. The mathematical framework of TO takes advantage of the metric invariance of Maxwell’s equations [3] to prove equivalence between change of coordinates and material properties [4]. Invisibility cloaks obtained analytically by TO require anisotropic, spatially varying permittivity and permeability. Such complex and extreme material parameters are not available in nature and have been achieved so far using metamaterials [5]. Experimental work have been carried out and proved the cloaking effect at microwave [6–8] and optical [9, 10] frequencies for both free space and carpet cloak configuration. However, the extreme material properties required for their fabrication are the main drawback for practical application of TO based cloaks, even if in some special cases one can relax the constraints on material parameters distribution [6] at the expense of an impedance mismatch. Moreover, most of the practical implementations rely on metamaterials to achieve such exotic material parameters through engineering their dispersion, which requires resonances and is linked to absorption through causality. The response of such cloaks is thus inherently narrowband and lossy, and the fabrication of such spatially varying metamaterials is a challenging task in the optical domain, which calls for alternative approaches.

Recently, a systematic methodology based on topology optimization have been proposed to realize all-dielectric cloaks with low index contrast [11], outlining a different strategy of design and fabrication. The cloaking effect is shown to appear at discrete angles of incidence for devices of increasing rotational symmetries. The same technique was applied for the design of graded index and binary cloaks working for both polarization [12]. From an experimental point of view, an optimized cloak made of single dielectric material has been demonstrated [13] in the microwave regime.

In this work we apply topology optimization to the design of cloaks made of Acrylonitrile butadiene styrene (ABS) dielectric material and concealing an ABS cylinder from a TE polarized line source excitation at microwave frequencies. This low index and low loss thermoplastic polymer is widely employed in 3D printing [14, 15], making the fabrication of the proposed devices straightforward. We assess the performance of various cloaks as a function of frequency by studying the wave scattering patterns and their reduction, as well as correlations of the total field with respect to the free space propagation case. Results from multi-objective optimization for one, two and three frequencies are presented and are shown to improve the bandwidth of operation. In this paper we propose a theoretical interpretation of the cloaking mechanism in terms of destructive interference between resonant modes of the system. This picture is validated numerically by quasimodal expansion method [16] and is shown to be radically different from the cloaking effect at stake in TO based structures.

2. Theory and numerical implementation

The problem stated in this paper is governed by Maxwell’s equation in time-harmonic regime (with convention exp(−iωt) discarded hereafter). Assuming that the parameters do not depend on z, the total electric field E_z excited by a line source located at r′= (x′, y′)^T obeys the following scalar wave equation:

\nabla \cdot (μ^{- 1} \nabla E_{z}) + k^{2} E_{z} = δ (r^{'}),

where ε and μ are the relative permittivity and permeability respectively, δ is the Dirac delta distribution, r = (x, y)^T is the spatial position vector and k = ω/c is the free space wavenumber. To solve this problem numerically, we use a Finite Element Method [17] (FEM) commercial package (COMSOL Multiphysics). A schema of the problem under study is depicted in Fig. 1. We put the source on the x axis and we model only half of the structure with homogeneous Neumann boundary condition on the symmetry plane (n·ΔE_z = 0). The free space is truncated by cylindrical Perfectly Matched Layers [18, 19] (PML) to damp propagating waves, with homogeneous Neumann condition on their outward boundaries.

Fig. 1 Numerical setup for FEM simulations (Comsol). Only one half of the domain is modelled with Perfect Magnetic Conductor (PMC) boundary conduction on the symmetry plane (TE polarization). The domain is truncated by Perfectly Matched Layers (PML).

Download Full Size | PDF

In the case of free space (ε = μ = 1), the solution of Eq. (1) is the Green’s function $G (r, r^{'}, k) = \frac{i}{4} H_{0}^{(1)} (k | r - r^{'} |)$ , where $H_{0}^{(1)}$ is the Hankel function of the first kind and $| r | = \sqrt{x^{2} + y^{2}}$ . The objective of the proposed cloak design is to minimize the diffracted field $E_{z d}^{c} = E_{z}^{c} - G$ in presence of the cloak and the scatterer in the background (air) domain Ω_B. For the case where the scatterer is not cloaked, the diffracted field is denoted $E_{z d}^{s} = E_{z}^{s} - G$ . Defining the integral quantities

W^{s} = \int_{Ω_{B}} {| E_{z d}^{s} |}^{2} d r, and W^{c} = \int_{Ω_{B}} {| E_{z d}^{c} |}^{2} d r,

the objective function to minimize is

ϕ_{1} (ε^{c}) = W^{c} / W^{s},

where ε^c is the permittivity in the design domain, subject to constraints ε_min ≤ ε^c ≤ ε_max.

To avoid small features and chessboard patterns, we put additional constraint on the gradient of permittivity. We define a second objective function

ϕ_{2} (ε^{c}) = \frac{h_{m}^{2}}{S} {| \nabla ε^{c} |}^{2},

where S = (R_c − R_s)²/2 is the area of the design domain and h_m is the maximum size of a mesh element.

Following the so-called solid isotropic material interpolation (SIMP) [20], we define a density function ξ such that

ε^{c} (ξ) = ε_{\min} + ξ^{p} (ε_{\max} - ε_{\min})

where p is a penalty factor (p = 1 here).

In order to enforce a binary design, we apply a threshold projection [21] at each iteration:

ξ (s) = \frac{\tanh (β ν) + \tanh (β (s - ν))}{\tanh (β ν) + \tanh (β (1 - ν))},

with ν = 0.5 and β = 2^N where N < 8 is the number of optimization iteration.

The topology optimization problem is thus to minimize the functional ϕ of the final design variable s

ϕ (s) = γ ϕ_{1} (s) + (1 - γ) ϕ_{2} (s),

subject to constraints 0 ≤ s ≤ 1 and where γ is a weight parameter between 0 and 1 (γ = 0.75 here). The gradient-based optimization code uses the Quasi-Newton method [22] to update the design iteratively. The sensitivity of the objective function with respect to s is calculated by the adjoint method [23]. The problem is discretized by the FEM using second order Lagrange elements and the solution is computed with a direct solver [24].

3. Results of optimization and study of cloaking performances

We optimized a permittivity distribution to cloak an ABS cylinder (ε_ABS = 2.67 − 0.0012i) of radius R_s = 3cm at a wavelength λ = 3cm from a line source located at x′= −2R_c = −12cm (cf. Fig. 1). The design domain is an annulus of thickness L_c = 3R_s = 6cm (R_s< r < R_s + L_c) and the permittivity is allowed to vary between ε_min = ε_air = 1 and ε_max = ε_ABS, starting from an initial empty configuration (s₀ = 0, i.e. $ε_{0}^{c} = 1$ ). The thickness of the annular background (air) domain and PML are L_b = R_s = 6 cm and L_PML = λ = 3cm respectively. Figure. 2(a) shows the optimized binary permittivity map obtained. The total electric field $E_{z}^{c}$ for the cloaked cylinder is shown on Fig. 2(c), and it can be seen that the scattering in free space is clearly reduced compared with the bare cylinder case (cf. Fig. 2(d)). This cloaking is confirmed by the plot of scattered field $E_{z d}^{c}$ (cf. Fig. 2(b)) which shows that the diffraction is extremely weak in the background domain and that the fields is concentrated within the cloak and the cylinder. Those field distributions indicates that there is a strong interaction between the cloaked object and the cloak, unlike TO-based cloaks where the wave is bent around the object, and do not penetrate in the cloaked region due to infinite material properties at the inner boundary of the cloak. This also means that the proposed optimized cloak is object dependant, but the same strategy can be applied to a PMC cylinder, so that any object hidden in it can be cloaked by the same cloak.

Fig. 2 (a) Optimized permittivity map ε^c(x, y). (b) Square norm of the scattered electric field for the cloaked case, showing very weak diffraction in free space outside the cloak. Real part of the total electric field with cloak (c) and for the bare cylinder (d).

Download Full Size | PDF

To confirm this qualitative observation, we computed the objective function ϕ₁ which quantifies the relative amount of energy scattered by the cloaked object compared with the uncloaked case. At the design wavelength (λ = 3cm), the scattering is reduced by 3 orders of magnitude as one can see on Fig. 3(a) (red solid line), but the operating range is quite narrow. Moreover, we have studied the influence of the position of the source on the performances and reported the results in Fig. 4. We found that the performances are maintained if the source is moved along the x axis on the left hand side of the cloak, but deteriorates rapidly as the source is off the x axis (the scattering reduction coefficient ϕ₁ increases by one order of magnitude for y′ ≃ λ/5, x′= −2R_c). In the limit where x′ → −∞, the excitation becomes a plane wave so we expect the cloak to work reasonably well in this case. This is confirmed by further computations in the case of a plane wave excitation coming from x < 0: the scattering reduction coefficient is ϕ₁ = 2.10⁻² at λ = 3cm (cf. Fig. 3(a), green dashed line), which is larger than the source point case but still represent a drastic reduction of 98% of the scattered field in comparison with a bare cylinder. More surprisingly, if the source is placed on the x-axis on the right hand side of the designed cloak, the scattering is still reduced significantly (cf. Fig. 4), resulting in a fairly good bi-directional cloaking.

Fig. 3 Study of cloaking performances. (a) Objective function ϕ₁ as a function of wavelength λ for source point (solid red line) and plane wave (green dashed line) excitation. (b) Correlation coefficient ρ as a function of wavelength λ for source point excitation for the real part (orange solid line) and imaginary part (dashed blue line) of the scattered electric field.

Download Full Size | PDF

Fig. 4 Scattering reduction coefficient ϕ₁ as a function of source position.

Download Full Size | PDF

Another parameter suitable to quantify the performances of the cloak is a correlation coefficient [25] calculated in the background between the electric fields for the cloaked and free space case, which is defined as:

ρ (X, Y) = \frac{E (X Y) - E (X) E (Y)}{\sqrt{E (X^{2}) - E {(X)}^{2}} \sqrt{E (Y^{2}) - E {(Y)}^{2}}}

where X is the Green’s function in free space G (real or imaginary part), Y is the electric field for the cloak and scatterer

E_{z}^{c}

(real or imaginary part) and

E (U) = \int_{Ω_{B}} U (r) d r

. We plotted on Fig. 3(c) the correlation coefficient as a function of wavelength, and one can observe that it is very close to unity at the design frequency for both real and imaginary part of the fields, which means that the electric field outside the cloak is close to the one propagating in free space, which indicates that the object is properly cloaked. Moving away from the design wavelength, ρ differs from unity which indicates that the performances deteriorates. This study of cloaking performances indicates that genuine topology optimization for a single frequency results in quite narrow bandwidth operation. To increase the frequency range of operation of the cloaking device, we propose to modify the optimization algorithm to several target frequencies.

4. Bandwidth enhancement with multi-frequency optimization

As an example, we set a number N_f of target frequencies equally distributed between λ = 2.7cm and λ = 3.3cm. The corresponding permittivity layouts are plotted in Fig. 5 (top) and it shows that as N_f increases, smaller features of dielectric map are required to reduce the scattering at higher frequencies. The resulting objective function for the different cases are reported in Fig. 5 (bottom) and it shows that the multi-frequency target leads to scattering around one order of magnitude higher than the single frequency case (red solid line). For N_f = 2, two dips in objective function appears, which shows the ability of the design process to optimize for multiple frequencies at the same time. For N_f = 3, the minimum scattering occurs at the centre of the frequency range and the bandwidth is clearly enhanced. The same conclusion applies for N_f = 5, even if the invisibility dip is slightly blueshifted. In this case, the bandwidth is 40% while it is of 33% bandwidth for the single frequency case. In summary, there is a tradeoff between bandwidth and cloaking performances. These results are analogous to those reported in [11], where optimization of cloaks were carried out with various incident plane wave angles in order to obtain omnidirectional devices. The optimized layouts with increased symmetry have shown better angular tolerance but at the expense of worse cloaking performances.

Fig. 5 Multifrequency optimization for bandwidth enhancement. (Top) Optimized permittivity distributions for N_f = 2, 3 and 5 frequencies (from left to right). (Bottom) Objective function ϕ₁ as a function of wavelength λ for various number of design frequencies N_f = 1 (red line), N_f = 2 (green line), N_f = 3 (orange line) and N_f = 5 (cyan line).

Download Full Size | PDF

5. Modal analysis and the origin of cloaking effect

Generally speaking, the diffractive properties of open systems such as the one we are dealing with in this paper can be studied at a more fundamental level by looking for both the generalized eigenfunctions and eigenvalues of the system: this is the so-called spectral problem. This boils down to solve for Eq. (1) without sources:

\nabla \cdot (μ^{- 1} \nabla ψ_{n}) + k_{n}^{2} ε ψ_{n} = 0.

Here ω_n = k_n c are the eigenfrequencies and ψ_n the quasimodes (also called leaky modes, quasi-normal modes or resonant state in the literature). Note that the eigenfrequencies are complex-valued because of the system being open, the real part being the resonant frequency and the imaginary part being related to dissipation processes (either absorption or leakage) and thus to the linewidth of the associated resonance.

The eigenfunctions are normalized in the following way:

〈 ε ψ_{n} | ψ_{n} 〉 : = \int_{Ω} ψ_{n} ψ_{m} d r .

A Quasimodal Expansion Method (QMEM) have been recently proposed as a tool for the analysis of photonic devices [16, 26]. The basic idea is to expand the total electric field $E_{z}^{c}$ solution of the problem with sources S onto the eigenmode basis:

E_{z}^{c} = \sum_{n} α_{n} ψ_{n} .

The coupling coefficients α_n indicates the relative contribution of each mode in the scattering process and are calculated as:

α_{n} = 〈 ε E_{z}^{c} | ψ_{n} 〉 = \frac{β_{n}}{k^{2} - k_{n}^{2}},

The coefficient β_n is given by the scalar product between the source and the eigenmode n as 〈β_n = S|ψ_n〉. Since in our case S = δ (r′) we have β_n = ψ_n(r′). One can then obtain a reduced order model by retaining M modes in the expansion:

E_{z}^{c} ≃ E_{z}^{c, rec} (M) = \sum_{n = 1}^{M} \frac{ψ_{n} (r^{'}) ψ_{n} (r)}{k^{2} - k_{n}^{2}} .

In this framework, the invisibility effect in such optimized cloak can be seen as a destructive interference between the contributions from the different resonant modes of the system, occurring for specific illumination conditions. The linear superposition of the different scattering channels cancels out in the background so that the diffracted field is dramatically reduced outside the cloak. We want to stress here that this mechanism is quite different from the cloaking produced by TO-based devices. Indeed, these ideal invisibility cloaks do not have any quasimode. The spectrum of such systems is composed of real eigenfrequencies associated with guided modes (which eigenfield is concentrated inside the cloaked object and cannot be excited by an external source), and of a continuous spectrum associated with radiation modes. A perfect TO cloak thus mimics effectively the spectrum of free space [27]. However practical realisation results in discretized cloaks based on metamaterials and leads to a system that exhibits leaky modes, so that the invisibility effect is narowband as well [28].

On the other hand, optimized dielectric cloaking system posses a much richer spectrum with an infinite countable number of leaky modes. The complex resonant interplay between these scattering channels is responsible for the cloaking effect at stake in that case. To illustrate this modal picture, we computed 1500 modes of the system and reconstructed the field using the expansion (11) with an increasing number of modes M. We define the reconstruction error as

e r r (M) = \int_{Ω} {| E_{z}^{c} - E^{c, rec} (M) |}^{2} d r,

and the reconstructed objective function as

ϕ_{1} (M) = \int_{Ω_{B}} {| E_{z d}^{c, rec} (M) |}^{2} d r .

where

E_{z d}^{c, rec} (M) = E_{z d}^{c, rec} (M) - G

is the reconstructed scattered field. These two quantities are plotted in Fig. 6 as a function of the number of modes M retained for the expansion. Firstly, one can see that the reconstruction error converges towards a low value (close to 10⁻⁴), which proves the validity of the QMEM. Secondly, we observe that the reconstructed objective function converges to the expected small value

ϕ_{1}^{opt}

as we add more and more modes to approximate the field, which illustrates the proposed modal picture. Note that this effect is similar to Fano resonances which stem from the interaction of a mode with a continuum of states. In our case the continuous spectrum is discretized by the use of PML and FEM but is fully included in the modal expansion (11): one can thus see the invisibility dip as an interaction between this continuous spectrum and all other leaky modes of the system.

Fig. 6 Convergence of the quasimodal expansion with the number of modes M. Reconstruction error err(M) (cyan solid line) and reconstructed objective function ϕ₁(M) (orange solid line).

Download Full Size | PDF

Finally, the cloaking here can also be analysed in the framework of scattering cancellation approach [29]: the role of the cloak is to cancel the contribution of multipole order, including higher order ones, even if the quasistatic approximation is not valid in our case. However, the practical implementation proposed here with structured dielectric differs from plasmonic cloaking that requires usually a thin cover of isotropic material with near-zero index [30].

6. Conclusion

We have reported the design of all-dielectric binary cloaking device using gradient based topology optimization, effectively reducing the scattering by a dielectric cylinder illuminated by a source point. A systematic analysis of has been provided to quantify the performance of the devices by using scattering reduction in comparison with the uncloaked case and correlation coefficient with respect to the free space propagation case. We extended the design to several frequencies and find that there is a trade-off between bandwidth and cloaking efficiency. Finally, we proposed a modal picture relying on destructive interferences to explain the invisibility effect. This study could pave the way to the design of practical cloaking devices and help the understanding of scattering reduction by complex shaped dielectric systems.

Acknowledgments

This work was funded by the Engineering and Physical Sciences Research Council (EPSRC), UK under a Programme Grant ( EP/I034548/1) “The Quest for Ultimate Electromagnetics using Spatial Transformations (QUEST)”.

References and links

1. J. B. Pendry, D. Schurig, and D. R. Smith, “Controlling electromagnetic fields,” Science 312, 1780–1782 (2006). [CrossRef] [PubMed]

2. U. Leonhardt, “Optical conformal mapping,” Science 312, 1777–1780 (2006). [CrossRef] [PubMed]

3. A. J. Ward and J. B. Pendry, “Refraction and geometry in Maxwell’s equations,” J. Mod. Opt. 43, 773–793 (1996). [CrossRef]

4. A. Nicolet, F. Zolla, Y. O. Agha, and S. Guenneau, “Geometrical transformations and equivalent materials in computational electromagnetism,” COMPEL 27, 806–819 (2008). [CrossRef]

5. H. Chen, C. T. Chan, and P. Sheng, “Transformation Optics and metamaterials,” Nat. Mater. 9, 387–396 (2010). [CrossRef] [PubMed]

6. D. Schurig, J. J. Mock, B. J. Justice, S. A. Cummer, J. B. Pendry, A. F. Starr, and D. R. Smith, “Metamaterial electromagnetic cloak at microwave frequencies,” Science 314, 977–980 (2006). [CrossRef] [PubMed]

7. B. Kanté, D. Germain, and A. de Lustrac, “Experimental demonstration of a nonmagnetic metamaterial cloak at microwave frequencies,” Phys. Rev. B 80, 201104 (2009). [CrossRef]

8. Y. Ma, Y. Liu, L. Lan, T. Wu, W. Jiang, C. K. Ong, and S. He, “First experimental demonstration of an isotropic electromagnetic cloak with strict conformal mapping,” Sci. Rep. 3, 2182 (2013). [CrossRef] [PubMed]

9. J. Valentine, J. Li, T. Zentgraf, G. Bartal, and X. Zhang, “An optical cloak made of dielectrics,” Nat. Mater. 8, 568–571 (2009). [CrossRef] [PubMed]

10. H. Chen and B. Zheng, “Broadband polygonal invisibility cloak for visible light,” Sci. Rep. 2, 255 (2012). [CrossRef] [PubMed]

11. J. Andkjær and O. Sigmund, “Topology optimized low-contrast all-dielectric optical cloak,” Appl. Phys. Lett. 98, 021112 (2011). [CrossRef]

12. J. Andkjær, N. Asger Mortensen, and O. Sigmund, “Towards all-dielectric, polarization-independent optical cloaks,” Appl. Phys. Lett. 100, 101106 (2012). [CrossRef]

13. L. Lan, F. Sun, Y. Liu, C. K. Ong, and Y. Ma, “Experimentally demonstrated a unidirectional electromagnetic cloak designed by topology optimization,” Appl. Phys. Lett. 103, 121113 (2013). [CrossRef]

14. B. Riddle, J. Baker-Jarvis, and J. Krupka, “Complex permittivity measurements of common plastics over variable temperatures,” IEEE Trans. Microw. Theory Techn. 51, 727–733 (2003). [CrossRef]

15. P. Deffenbaugh, R. Rumpf, and K. Church, “Broadband microwave frequency characterization of 3-d printed materials,” EEE Trans. Compon. Packag. Manuf. Techno. 3, 2147–2155 (2013). [CrossRef]

16. B. Vial, F. Zolla, A. Nicolet, and M. Commandré, “Quasimodal expansion of electromagnetic fields in open two-dimensional structures,” Phys. Rev. A 89, 023829 (2014). [CrossRef]

17. J. Jin, The Finite Element Method in Electromagnetics (John Wiley & Sons, 2014).

18. J.-P. Berenger, “A perfectly matched layer for the absorption of electromagnetic waves,” J. Comput. Phys. 114, 185–200 (1994). [CrossRef]

19. M. Lassas, J. Liukkonen, and E. Somersalo, “Complex Riemannian metric and absorbing boundary conditions,” J. Math. Pure. Appl. 80, 739–768 (2001). [CrossRef]

20. M. P. Bendsøe and O. Sigmund, “Material interpolation schemes in topology optimization,” Archive of Applied Mechanics 69, 635–654 (1999). [CrossRef]

21. F. Wang, B. Lazarov, and O. Sigmund, “On projection methods, convergence and robust formulations in topology optimization,” Struct. Multidiscip. O. 43, 767–784 (2011). [CrossRef]

22. R. Tapia, “Diagonalized multiplier methods and quasinewton methods for constrained optimization,” Journal of Optimization Theory and Applications 22, 135–194 (1977). [CrossRef]

23. J. S. Jensen and O. Sigmund, “Topology optimization of photonic crystal structures: a high-bandwidth low-loss t-junction waveguide,” J. Opt. Soc. Am. B 22, 1191–1198 (2005). [CrossRef]

24. O. Schenk, A. Wachter, and M. Hagemann, “Matching-based preprocessing algorithms to the solution of saddle-point problems in large-scale nonconvex interior-point optimization,” Computational Optimization and Applications 36, 321–341 (2007). [CrossRef]

25. D. Bao, R. Mitchell-Thomas, K. Rajab, and Y. Hao, “Quantitative study of two experimental demonstrations of a carpet cloak,” IEEE Antennas Wireless Propag. Lett. 12, 206–209 (2013). [CrossRef]

26. C. Sauvan, J. P. Hugonin, I. S. Maksymov, and P. Lalanne, “Theory of the spontaneous optical emission of nanosize photonic and plasmon resonators,” Phys. Rev. Lett. 110, 237401 (2013). [CrossRef] [PubMed]

27. M. Farhat, P.-Y. Chen, S. Guenneau, S. Enoch, R. McPhedran, C. Rockstuhl, and F. Lederer, “Understanding the functionality of an array of invisibility cloaks,” Phys. Rev. B 84, 235105 (2011). [CrossRef]

28. A. Rahmani, M. J. Steel, and P. C. Chaumet, “Invisibility and supervisibility: Radiation dynamics in a discrete electromagnetic cloak,” Phys. Rev. B 87, 045430 (2013). [CrossRef]

29. P.-Y. Chen, J. Soric, and A. Alù, “Invisibility and cloaking based on scattering cancellation,” Advanced Materials 24, OP281–OP304 (2012). [PubMed]

30. A. Alù and N. Engheta, “Achieving transparency with plasmonic and metamaterial coatings,” Phys. Rev. E 72, 016623 (2005). [CrossRef]

Topology optimized all-dielectric cloak: design, performances and modal picture of the invisibility effect

Abstract

1. Introduction

2. Theory and numerical implementation

3. Results of optimization and study of cloaking performances

4. Bandwidth enhancement with multi-frequency optimization

5. Modal analysis and the origin of cloaking effect

6. Conclusion

Acknowledgments

References and links

Cited By

Figures (6)

Equations (15)

Optics Express