## Abstract

We propose and develop a technique for designing a special class of nonmagnetic metamaterials possessing desired dielectric and optical properties over a broad frequency band. The technique involves the design of nanostructured metallodielectric materials (photonic crystals) with a special layered geometry where the metal content in each layer has to be determined using a fitting procedure. For illustration, we demonstrate the performance of our technique for tailoring metamaterials having epsilon-near-zero and on-demand refractive index (real or imaginary part) over a frequency band. One-, two-, as well as three-dimensional geometries have been considered. In the one-dimensional and two-dimensional cases, the results of semi-analytical calculations are validated by *ab initio* FDTD simulations.

©2013 Optical Society of America

## 1. Introduction

In the last decade, the design and fabrication of various man-made materials with bizarre optical properties not available in nature became a very promising field of research and industry. Spectacular progress in this field has been achieved due to an unprecedented freedom in controlling the inner structure of materials on the micro- and nano-scale resulted from advances in modern fabrication technology. Nevertheless, it seems that in many cases current materials science is not able to provide the industry with technology solutions on how to tailor metamaterials (MMs) possessing needed mechanical, elastic, thermal, electrical, and optical properties. In other words, practical demands outrun our production capabilities. As a result, a lot of projects remain unrealized. One glowing example is the project of “invisibility” cloak, especially interesting for broadband optical applications, whose realization is closely allied to a new metamaterial design [1].

In this paper, by combining materials with the positive permittivity (dielectric) and negative permittivity (metal), we develop an easy-to-use technique for the design of nonmagnetic MMs having desired optical properties over a frequency band that in turn may facilitate the design of new nano-devices and nano-components. In essence, we propose a recipe for designing a special class of MMs: broadband MMs, which allow one to extend the range of possible applications which could play an important role in various branches of science and engineering. Recently, novel MMs exhibiting broad bandwidth and low losses in both the microwave [2] and optical [3] range have been proposed and fabricated. Nevertheless, one can state that such MMs are still in the early stage of their development. We foresee that many fields such as communications, optical processing, data storage, microscopy, medicine, etc. could ultimately benefit from the development of these MMs.

As a starting point, we consider a stratified metallodielectric composite. More specifically, we deal with a multilayered composite translationally invariant in the *x*-direction (the direction of propagation). In the *z*-direction, the composite consists of identical unit cells of the width *S _{z}*. Taken together, the cells form a one-dimensional (1d) photonic crystal. In turn, each unit cell may be considered as a gradient film, i.e., allowance is made for the varying permittivity

*ε*(

*z*) within each cell in the

*z*-direction.

A general consideration of the optical properties of the gradient film can be given in terms of the wave equation. However, a closed-form analytical solution of this equation, being dependent on the specific form of the profile *ε*(*z*), is possible only in special cases [4]. As to the solution of the inverse problem (determination of the function *ε*(*z*) from the predetermined optical properties), which is the primary focus of our work, it presents a challenge and calls for even more sophisticated technique. Due to that, here we suggest another approach.

In some cases, as in our previous work on broadband epsilon-near-zero (ENZ) MMs [5], the desired effect can be achieved by means of the direct specification of the permittivity profile. That work has been based on the fact that for a layered medium, if the electric field is normal to the layers and each layer has its local permittivity close to zero at a certain frequency, the effective permittivity of the whole medium at this frequency must be close to zero, too. Nulling the permittivities of particular layers at frequencies within a band, we hence can null the effective permittivity over this band. However, such an approach, as well as a somewhat modified approach [6–8], is inapplicable if we want to obtain the permittivity different from zero or on-demand refractive index over a finite frequency band. The present approach is technically simpler, because, as will be seen, it involves no thickness distribution function, and more general, because it allows one to design broadband MMs with various optical properties [9,10]. Here, we consider it in detail addressing such issues as its accuracy and some possible extensions to the three-dimensional case.

The approximation used in our work is in essence the conventional mean-field approach because it does not take into account such effects as magnetoelectric coupling and spatial dispersion. The role of the above effects is currently undergoing revision (see, e.g., Refs [11,12]. and references therein). Here, these effects are considered to be weak that usually takes place for small wave vectors and off-resonance regime that, as will be shown, is our case. At the same time, the accuracy of our approximations for some specific cases has been estimated with the use of *ab initio* finite-difference time-domain (FDTD) simulation and it is sufficient for many potential applications. More accurate technique taking into account local-field effects will be considered elsewhere.

To conclude this section, it must be emphasized that our approach is completely analytical. This can make it advantageous as compared to usually more complicated numerical techniques used for description of the optical properties of metallodielectric photonic crystals [13–15] which often involve problems associated with the poor convergence and a prohibitive computational cost.

## 2. Geometry and technique

First, we assume that the electric field is applied along the *z* axis. Furthermore, we deal with the long-wavelength regime in which *S _{z}*<<

*λ*/Re(

*n*) where

_{eff}*n*is the

_{eff}*effective*refractive index of the composite for the chosen orientation of the electric field and

*λ*is the vacuum wavelength of the light. So, we do not take into consideration the effect of the finite size of the unit cells (retardation effects) which can result in the non-zero phase advance across the cell [2]. Then, the effective permittivity of the composite

*ε*= (

_{eff}≡${\epsilon}_{eff}^{z}$*n*)

_{eff}^{2}may be written as

*ε*(

*z*) may be approximated as a step-like (piece-wise) function. In this case each unit cell is considered to be consisting of

*N*thin parallel layers (see, e.g., Fig. 1 from [5]), and the effective permittivity may be written as the harmonic averageEquation (1a) may be derived from Eq. (1b) in the limit of the very thin layers or, equivalently, as N→∞. Equations (1a) and (1b) are frequently used as basic equations for designing macroscopically homogeneous materials with needed useful properties [16].

In fact, there are many ways to obtain the desired profile *ε*(*z*). For example, when dealing with metallodielectric composites, both the metal plasma frequency and collision frequency can be tuned by imposing a gradient of temperature in the *z*-direction that results in the corresponding alteration of the permittivity profile along *z* [17]. Another way is to deal with a suspension of small metallic spheres in a dielectric host when the particle volume fraction is a function of *z* [18] or, alternatively, with a suspension of small dielectric particles in a metal host [19]. Intermetallics and metal alloys can be also used for designing particular materials with desired optical properties [20]. A good list of literature on how to practically realize designed profiles *ε*(*z*) is given in the review by Shvartsburg et *al* [4]. There is, however, a conceptually different approach which deserves attention to design nanostructured MMs possessing desired broadband dielectric and optical properties.

Strictly speaking, Eq. (1) is valid only in one-dimensional (1d) case and within the framework of the formalism of the local permittivity. However, as an approximation, it can be used also for more complicated systems (two-phase metallodielectric composites) with properly chosen geometry. The key idea is that this geometry can be chosen in such a way that the applied electric field remains *almost* parallel to the metal/dielectric interface within the unit cell width. To explain how to realize such a geometry let us look at Fig. 1. Here, we sketch two metal/dielectric unit cells; each cell may be considered as consisting of 8 parallel layers of the uniform thickness. In its turn, each layer may be represented as a set of parallel segments oriented along the *z* axis. Along the *y* axis, the structure is considered to be homogeneous. The segments in the form of dielectric channels (pores) in a metallic host, as well as in the form of metallic channels (wires) in a dielectric host are allowed. In both cases, the width of the channels may differ in different layers. We note that all layers are considered to be of uniform thickness in contrast to our previous work [5], where we have introduced a thickness distribution function.

Because the local permittivity changes now along the *x*-axis as well, the composite becomes two-dimensional (2d) rather than 1d. However, if the size of the unit cell in the *x-*direction (*S*_{x}) is small (*S*_{x}<<*S*_{z}), the composite may be considered as quasi one-dimensional. Furthermore, due to the continuity of the electric field which is tangential to the metal/dielectric interfaces, the permittivity of the i*th* layer may be written as

*ε*

_{1}and

*ε*

_{2}are the permittivities of the metallic and dielectric phases,

*f*is the volume fraction of the metal phase in the i

_{i}*th*layer, and

*f(z)*is the corresponding distribution function of the metal phase. Substitution of Eq. (2a) into Eq. (1b) and Eq. (2b) into Eq. (1a) gives

*s´*=

*z*/

*S*, and

_{z}*s*=

*ε*/(

_{2}*ε*–

_{1}*ε*).

_{2}Looking at Eqs. (3) for the effective permittivity, it is easy to see that it has no poles. In other words, within the framework of our consideration, this component of the effective permittivity is nonresonant. At the same time, it can take zero values; when neglecting the imaginary part of the phase permittivities, zeros of the function *ϵ _{eff}* occur at

*f*(

*z*) = -

*s*(or

*f*= -

_{i}*s*). This feature of the effective permittivity results from a specific geometry of the photonic crystals under consideration where localized surface plasmons cannot be excited and hence resonances are not present.

Thus, to tailor a metamaterial with the desired optical properties, the problem may be reduced to minimizing the functional of the form

*F*is an objective function, and the sought function

*f*has to be determined using a proper minimization procedure. Thus, in the present approach the unknown function is the distribution of metal content along z-coordinate, while in our previous work [5] we look for the distribution function in terms of the layer thickness.

## 3. Results and discussion

To demonstrate the potential of our technique, let us consider some examples. In this section, we have set *ε*_{2} = 2.5 for the permittivity of a dielectric and used the Drude representation of the form ${\epsilon}_{1}\left(\omega \right)=1-1/\left[{\omega}^{2}/{\omega}_{p}^{2}+i\omega \gamma /{\omega}_{p}^{2}\right]$ with $\gamma /{\omega}_{p}=0.01$ for the permittivity of a metal. In our examples we have found metal content in each layer by optimizing the optical characteristics of interest over the frequency bands of $\left(0.32-0.4\right){\omega}_{p}$ and $\left(0.35-0.4\right){\omega}_{p}$. For these frequency bands, the metal permittivity changes from −8.76 at $\omega /{\omega}_{p}=0.32$ and −7.2 at $\omega /{\omega}_{p}=0.35$ to −5.2 at $\omega /{\omega}_{p}=0.4$that corresponds to the visible range for such metals as gold and silver.

#### 3.1 Design of epsilon-near-zero MMs

ENZ MMs have been in the focus of scientific interest for the last few years due to their versatile potential applications [21] which include those using cloaking effect [22]. However, so far the design of *broadband* ENZ MMs remains a challenging task. A possible solution has been proposed in our previous work [5]. Our present approach may be considered as an alternative.

The objective function for this case may be represented as $F={\left[\mathrm{Re}{\epsilon}_{eff}\left(\omega ,f\right)\right]}^{2}$. For the geometry shown in Fig. 1, a near-zero *z*-component of the dielectric tensor is obtained after dividing each unit cell into 8 layers and finding the metal content in each of them. Then, using Eq. (3a), the effective permittivity (see Fig. 2
) is calculated. Close inspection shows that |Re*ε _{eff}*| does not exceed 0.05 in the actual frequency band of $\omega /{\omega}_{p}$ = 0.35−0.4. As shown in Fig. 2, Re

*ε*is a function oscillating about zero within the actual band where the number of its zeros corresponds to the number of layers in the unit cell. In this example, our model for Re

_{eff}*ε*yields 7 zeros within the band. It is so because in this particular case the metal contents in two of eight layers are very close. To assess the accuracy of fitting, we have also determined, and shown in Fig. 2, the root mean square of the deviation of the fit which in this case may be defined as $rms=\sqrt{\left[{\displaystyle \sum _{i=1}^{M}\mathrm{Re}{\epsilon}_{eff}\left({\omega}_{i}\right)}\right]/M}$with

_{eff}*ω*running through the band [

_{i}*ω*

_{1},

*ω*

_{2}].

#### 3.2 Design of MMs with on-demand refractive index

Among MMs with on-demand refractive index, of particular interest are those with low and ultra-low refractive index (ULRI) MMs. They are important for many microelectronic [23] and optical [24–26] applications.The objective function may be represented as $F={\left[\mathrm{Re}\sqrt{{\epsilon}_{eff}\left(\omega ,f\right)}-{n}_{d}\right]}^{2}$. We have considered the effective refractive index of a metallodielectric composite whose Re*n _{eff}* has been fitted to

*n*= 0.5. For the case when each unit cell is divided into 16 layers, our results for two frequency bands are shown in Fig. 3 . As in the previous case, Re

_{d}*n*oscillates within the actual band. After broadening the band, the oscillation amplitude builds up, especially at low frequencies. As can be seen, at

_{eff}*N*= 16 the accuracy of fitting reaches about 1% (i.e., |Re

*n*-0.5|<0.005) over the band $\left(0.35-0.4\right){\omega}_{p}$. After broadening the frequency band, the accuracy becomes worse and reaches about 7% on the low-frequency side of the band. It should be noted that a lower bound exists on the allowable value of

_{eff}*n*; in other words, it cannot be made arbitrary close to zero. Indeed, $\mathrm{Re}{n}_{eff}=\sqrt{\left[\mathrm{Re}{\epsilon}_{eff}+\sqrt{{\left(\mathrm{Re}{\epsilon}_{eff}\right)}^{2}+{\left(\mathrm{Im}{\epsilon}_{eff}\right)}^{2}}\right]/2}$. Because both Re

_{d}*ε*and Im

_{eff}*ε*are limited within the considered range, Re

_{eff}*n*is limited as well: $\mathrm{Re}{n}_{eff}\ge \sqrt{\left[\mathrm{Re}{\epsilon}_{eff}^{\mathrm{min}}+\sqrt{{\left(\mathrm{Re}{\epsilon}_{eff}^{\mathrm{min}}\right)}^{2}+{\left(\mathrm{Im}{\epsilon}_{eff}^{\mathrm{min}}\right)}^{2}}\right]/2}$. Thus, for obtaining a small refractive index, the imaginary part of the effective permittivity must be as small as possible while its real part must be large and negative.

_{eff}#### 3.3 Design of MMs with high absorption efficiency

Highly absorbing (blackbody-like) optical materials have received considerable attention in recent years, especially due to their possible use for various optical instruments and sensors, coherent thermal emitters, photodetectors, photovoltaic devices, *etc* [27–30].

To design a metamaterial with the constant absorption coefficient *α _{c}* or imaginary part of the refractive index

*k*over a frequency band, the objective function may be represented as $F={\left[\frac{\omega}{c}\mathrm{Im}\sqrt{{\epsilon}_{eff}\left(\omega ,f\right)}-{\alpha}_{c}\right]}^{2}$or $F={\left[\mathrm{Im}\sqrt{{\epsilon}_{eff}\left(\omega ,f\right)}-{\kappa}_{c}\right]}^{2}$, respectively. In Fig. 4 we show the imaginary part of the refractive index fitted to

_{c}*k*= 0.45 for two frequency bands, $\left(0.35-0.4\right){\omega}_{p}$and $\left(0.32-0.4\right){\omega}_{p}$, where each unit cell is assumed to consist of 16 layers. In contrast to the case of Re

_{c}*n*, oscillations of Im

_{eff}*n*occur on the high-frequency side of the band.

_{eff}For completeness, in Fig. 5
we give the distribution function of the metal phase *vs* the vertical coordinate for three above cases, Re*ε _{eff}* = 0 (8 layers within the unit cell), Re

*n*= 0.5 (16 layers within the unit cell) and Im

_{eff}*n*= 0.45 (16 layers within the unit cell).

_{eff}## 4. Accuracy of the approximation used

A natural question arises regarding the accuracy of the proposed method. In fact, Eq. (1) is exact in the 1d case. As to Eq. (2), it is a remarkable fact that for metal alloys their effective permittivity can be often considered as a weighted average of the permittivities of the parent metals (see, e.g., [31,32] and references therein); in general, it is shown to be so for thermodynamically ideal binary mixtures [33]. To check the accuracy of the approximation used in the 2d case for the geometry shown in Fig. 1, we have performed *ab initio* numerical simulation using FDTD technique and compared the results with those obtained with the use of Eq. (3a). The simulation for the 1d case has been also performed for the comparison. To be more specific, we have found the effective permittivity with the use of the computed *z*-component of the local electric field,

*V*is the unit cell volume, and then compared it with that obtained with the use of Eq. (3a) for the above ENZ MM. The number of the unit cells in the

*x*direction is 16, the sizes of the unit cell are

*S*≈0.0094 λ and

_{x}*S*≈0.066 λ with respect to the band center. It should be noted that such a choice of the unit cell size satisfies the typical evaluation criteria for the correct homogenization of MM lattices [34],(it should be reminded that we deal with off-resonant regime). In other words, the effects of spatial dispersion (that, generally, can be important for photonic crystals [12]) are considered to be weak. Besides, to obtain the physically correct time-independent effective permittivity, the necessary conditions for the steady state regime have been provided.

_{z}Our results for the effective permittivity *ε _{eff}* are given in Fig. 6
. As can be seen, the agreement between our approximation and FDTD simulation is very good in the 1d case, but with somewhat smoothed out oscillations of Re

*ε*curve. In the 2d case, the agreement is fairly good: it is perfect at high frequencies and certainly worse at low frequencies. Besides,some oscillations are not resolved in the actual range. Within the band $\left(0.35-0.4\right){\omega}_{p}$ the deviation of Re

_{eff}*ε*from zero is about 0.1. The agreement is much better for the 1d geometry than in the 2d case. Obviously, this happens because Eq. (1) is exact (for the 1d case), while Eq. (3) is approximate (for the 2d case).

_{eff}It should be noted that FDTD simulations never reproduce the results of rigorous calculations exactly, but only to within a simulation accuracy. We believe that the disagreement between our approach and FDTD simulation in the 1d case can be attributed mainly to the inaccuracy of the simulation procedure used and to some extent to the spatial dispersion. The agreement can be improved by decreasing the Yee cell size that specifies the level of discretization (in our case, it is about 7.5 × 10^{−5} λ) and by increasing the thickness of the perfectly matched layers (we have taken that about one wavelength corresponding to the band center) and the total film thickness. However, as we have made sure, in the 2d case the change of the computational parameters (which can result in a significant increase in the computation time) has only a moderate effect on the agreement. Furthermore, it is known that when dealing with ultrathin films, the effective permittivity for the film as a whole may differ from that for the corresponding bulk (infinite) medium [35]. This results from the surface effect (due to different environment of inner and surface layers, averaged (effective) fields in these layers are also different). Generally, the number of layers, which are necessary to achieve the bulk regime, should depend on the coupling between adjacent layers [36]. As was shown for coupled-fishnet structures, when the distance between adjacent layers is small and hence the coupling is strong, at least four layers are required [37]. The surface effect is pronounced in the case of a few layers only while we deal with 16 unit cells along the *x*-direction. Besides, due to a weak coupling between the cells, it should not be strong as the electric field is oriented mainly along the boundaries (along the wires). The effect of the finite number of layers of 1d and 2d lattices with periodic boundary conditions in the directions normal to the propagation direction on the effective permittivity has been considered using the transfer matrix method [38]. Dealing with a wire medium, the authors have found that only if it consists of one layer, a noticeable deviation of the obtained effective permittivity from the bulk one occurs, while for the three and five layers the deviation is negligible. The properties of the surface plasmon polaritons of a layered metallodielectric MM with the use of the effective medium theory have been examined in [39]. In this case, one should consider two diagonal components of the effective permittivity tensor. The authors have revealed a close proximity between their results for the two-layer medium and the bulk case. So we believe that the main reason for the above disagreement in the 2d case is the inaccuracy of theapproximation used. As was noted, Eq. (1) is valid for the 1d case. In fact, as we move from one layer to adjacent parallel layer, the electric field (which, we remind, is polarized along the *z* axis) undergoes a jump discontinuity, i.e., strictly speaking, is not continuous. This discontinuity could be reduced by increasing the layer thickness (however, it is necessary to keep in mind that this thickness cannot be too large because the total thickness of the unit cell should be much less than wavelength) or by decreasing the size of the wires or pores in the lateral direction, i.e., *S*_{x}*/S*_{z} ratio. Mathematically, as can be seen from Fig. 1, a necessary condition of parallelism between the electric field and the metal/dielectric interfaces may be written in the form

*ε*| becomes very small [40].

_{eff}In the above calculations, we have dealt with the relatively narrow frequency bands $\omega /{\omega}_{p}=0.35-0.4$ and $\omega /{\omega}_{p}=0.32-0.4$. However, good results may be also obtained for broader bands. To do this, we should consider the unit cell with a larger number of layers.

In many cases, small values of the imaginary part of the effective permittivity or refractive index are important in addition to the desired values of their real part over an actual frequency band. While neither Im*ε _{eff}* nor Im

*n*can be nulled, there are at least two ways to reduce them dealing with two-phase metamaterials. One way is to use metal having smaller values of Im

_{eff}*ε*

_{1}; for example, for silver below 3 eV, where the permittivity follows the Drude model, the dimensionless damping parameter may be estimated as $\gamma /{\omega}_{p}\approx 0.01$. Recently, the alkali-noble intermetallics have been shown to hold promise as low-loss plasmonic materials [41]. Another way to reduce Im

*ε*is to use a dielectric phase with smaller value of its permittivity.

_{eff}The results reported here were obtained for the fixed value of the permittivity or refractive index. It is, however, clear that our approach can be used for design of graded-index (or graded-permittivity) films and coatings. So, when the profile of the refractive index matches the refractive index of air, on the one side, and the refractive index of substrate on the other side, broadband antireflective coatings (being of interest, in particular, for improving the efficiency of solar cells) can be realized. Furthermore, graded-index materials having anywhere the refractive index smaller than unity play the critical role for realization of optical cloaking [26].

## 5. Generalization to the three-dimensional case

In the 3d case, there are many options for designing MMs possessing the needed optical properties within a frequency band. Direct approach dealing with two-phase composites is to use the widely known Bergman-Milton spectral representation for the effective permittivity [42]. Such representation has been recently applied for designing broadband ENZ MMs with the use of a simple mixing rule for the effective permittivity as an average of the permittivities of constituents weighted by their volume fractions [6]. It should be, however, noted that, similar to our previous consideration, the above results are feasible to design only one of the diagonal components of the effective permittivity tensor. In other words, designed MM is anisotropic, and while one of the components of the tensor of the effective permittivity is properly fitted, two other components can be arbitrary. This results from a specific geometry of MM (a flat metal-dielectric layered stack). At the same time, consideration can be given to the 3d case when any of the diagonal components of the tensor can be properly tailored.

According to Bergman and Milton, the effective permittivity of a two-phase three-dimensional isotropic composite may be represented as

*G*(

*x*), containing all information about geometry of the composite, satisfies the sum rules ${\int}_{0}^{1}G(t)dt}=1$ and ${\int}_{0}^{1}tG(t)dx}={p}_{2}/3$, and

*p*

_{1}and

*p*

_{2}are the volume fractions of the phases 1 and 2, respectively.

With a knowledge of the effective permittivity within a broad frequency range, we could extract the spectral density, solving an inverse problem [43–46]. However, a problem remains how to relate the function *G*(*x*) to the real geometry of the composite. This problem is not trivial, and the above relationship is known in some particular cases only (see, e.g., [47]); in general, geometrical information is shown to be incorporated into the spectral density function via its moments [48]. Here we do not consider this problem in more detail and redirect our attention to other possibility specifying the internal geometry of the composite.

To determine the effective permittivity, the Maxwell-Garnett approximation (MGA) is shown to be rather accurate for rarefied suspensions of spherical particles (see, e.g., [49–51] and references therein). For the case of multiphase suspensions, MGA is of the form [52]

*p*is the volume fractions of spherical particles of i

_{i}*th*kind (in the following all

*p*are taken to be equal to reduce the number of free parameters). Each particle is considered to be an alloy of two metals with the permittivity ${\epsilon}_{1}\left(\omega \right)$ and ${\epsilon}_{2}\left(\omega \right)$. For simplicity, we assume that both ${\epsilon}_{1}\left(\omega \right)$ and ${\epsilon}_{2}\left(\omega \right)$are of the Drude form but with different plasma frequencies ${\omega}_{p1}$ and ${\omega}_{p2}$such that ${\omega}_{p1}^{2}/{\omega}_{p2}^{2}=0.8/1.2$:${\epsilon}_{1}\left(\omega \right)=1-0.8/\left[{\omega}^{2}/{\omega}_{p*}^{2}+i\omega \gamma /{\omega}_{p*}^{2}\right]$ and ${\epsilon}_{2}\left(\omega \right)=1-1.2/\left[{\omega}^{2}/{\omega}_{p*}^{2}+i\omega \gamma /{\omega}_{p*}^{2}\right]$ where ${\omega}_{p*}^{2}=\left({\omega}_{p1}^{2}+{\omega}_{p2}^{2}\right)/2$. In its turn, the effective permittivity of the binary alloy particles is considered to be of the form (2), and the permittivity of a dielectric host ${\epsilon}_{0}=3.1$. Thus, the problem can be reduced to minimization of the functional (4) with the effective permittivity specified by Eq. (9), where the sought parameters are the volume fractions of the metal phase 1 (or 2) in the alloy particles.

_{i}Figure 7
shows the results of our computations for $\mathrm{Re}{\epsilon}_{eff}\left(\omega \right)$fitted to zero over the frequency band $\left(0.35-0.4\right){\omega}_{p*}$ for a rarefied suspension of metal alloy particles of 10 kinds for 3 values of the damping parameter. The total volume fraction of the particles is 0.06 (6%). As can be seen, the obtained effective permittivity oscillates around zero over the band, and an increase in the damping parameter allows one to sufficiently decrease the amplitude of oscillations and hence to reduce *rms*. In particular, after increasing the damping parameter by the factor 1.6, the root mean square of the deviation from zero becomes five times smaller.Because *ε* = 0 is the condition for plasmon excitation, our results show that bulk plasmons can be excited in a frequency range, as well as surface plasmons. Thus, longitudinal electromagnetic waves can propagate in such broadband MMs within a band.

## 6. Concluding remarks

We have shown that by specifying a gradient of the local permittivity of a MM along the direction of the applied electric field, it is possible to design a structure of the unit cell of 1d or 2d photonic crystal in such a way that the crystal effective permittivity or refractive index would have a predetermined value over a frequency band. To do that, we derive an analytical representation for the effective permittivity of MMs of a special kind, namely, nonresonant subwavelength photonic crystals. In addition, by proper choice of the objective function we can also specify the needed profile of the effective permittivity or refractive index in the band. For illustration, we have shown how our approach works for designing broadband epsilon-near-zero, on-demand refractive index, and blackbody-like MMs. The validity of our approach has been supported with the use of *ab-initio* FDTD simulations and some limitations on its applicability have been discussed. In particular, the approach is shown to be accurate in the 1d case and fairly good in the 2d case. Finally, we have shown how our technique can be extended to the 3d case dealing with rarefied suspensions of metal alloy spheres.

One of the issues impeding the further development of plasmonic MMs generally and considered broadband MMs specifically is their relatively high dielectric losses. For some applications, however, high losses can be desired [53]. At the same time, there is, in principle, an opportunity to significantly reduce dielectric losses of broadband MMs using loss compensation by gain media [7,54] or superconducting constituent materials [55].

At present, the fabrication of broadband MMs appears to be a challenging task. Nevertheless, taking into account the rapid progress in nanotechnology, we are hopeful that this challenge will be overcome in the near future, first of all in the IR and microwave range.

## Acknowledgments

The authors thank National Center for High-Performance Computing for computer time and facilities. The work was supported by NCKU Project of Promoting Academic Excellence and Developing World Class Research Centers. V.U.N. acknowledges partial support from National Science Council, Taiwan, Grant No. 100-2112-M-001-025-MY3.

## References and links

**1. **W. Cai and V. Shalaev, *Optical Metamaterials: Fundamentals and Applications* (Springer, 2010).

**2. **R. Liu, Q. Cheng, J. Y. Chin, J. J. Mock, T. J. Cui, and D. R. Smith, “Broadband gradient index microwave quasi-optical elements based on non-resonant metamaterials,” Opt. Express **17**(23), 21030–21041 (2009). [CrossRef] [PubMed]

**3. **J. Valentine, S. Zhang, T. Zentgraf, and X. Zhang, “Development of bulk optical negative index fishnet metamaterials: achieving a low loss and broadband response through coupling,” Proc. IEEE **99**(10), 1682–1690 (2011). [CrossRef]

**4. **A. V. Shvartsburg, V. Kuzmiak, and G. Petite, “Optics of subwavelength gradient nanofilms,” Phys. Rep. **452**(2-3), 33–88 (2007). [CrossRef]

**5. **A. V. Goncharenko and K. R. Chen, “Strategy for designing epsilon-near-zero nanostructured metamaterials over a frequency range,” J. Nanophotonics **4**(1), 041530 (2010). [CrossRef]

**6. **L. Sun and K. W. Yu, “Strategy for designing broadband epsilon-near-zero metamaterials,” J. Opt. Soc. Am. B **29**(5), 984–989 (2012). [CrossRef]

**7. **L. Sun and K. W. Yu, “Strategy for designing broadband epsilon-near-zero metamaterial with loss compensation by gain media,” Appl. Phys. Lett. **100**(26), 261903 (2012). [CrossRef]

**8. **L. Sun, K. W. Yu, and X. Yang, “Integrated optical devices based on broadband epsilon-near-zero meta-atoms,” Opt. Lett. **37**(15), 3096–3098 (2012). [CrossRef] [PubMed]

**9. **A. V. Goncharenko, V. U. Nazarov, and K. R. Chen, “Metallodielectric broadband metamaterials,” SPIE Newsroom (Feb. 6, 2012). DOI: 10.1117/2.1201201.0040207. [CrossRef]

**10. **A. V. Goncharenko, V. U. Nazarov, and K. R. Chen, “Development of metamaterials with desired broadband optical properties,” Appl. Phys. Lett. **101**(7), 071907 (2012). [CrossRef]

**11. **A. Alù, “First-principles homogenization theory for periodic metamaterials,” Phys. Rev. B **84**(7), 075153 (2011). [CrossRef]

**12. **C. R. Simovski, “On electromagnetic characterization and homogenization of nanostructured metamaterials,” J. Opt. **13**(1), 013001 (2011). [CrossRef]

**13. **K. Sakoda, N. Kawai, T. Ito, A. Chutinan, S. Noda, T. Mitsuyu, and K. Hirao, “Photonic bands of metallic systems. I. Principle of calculation and accuracy,” Phys. Rev. B **64**(4), 045116 (2001). [CrossRef]

**14. **Z. Y. Li and K. M. Ho, “Analytic modal solution to light propagation through layer-by-layer metallic photonic crystals,” Phys. Rev. B **67**(16), 165104 (2003). [CrossRef]

**15. **A. A. Krokhin, E. Reyes, and L. Gumen, “Low-frequency index of refraction for a two-dimensional metallodielectric photonic crystal,” Phys. Rev. B **75**(4), 045131 (2007). [CrossRef]

**16. **V. U. Nazarov, “Bulk and surface dielectric response of a superlattice with an arbitrary varying dielectric function: A general analytical solution in the local theory in the long-wave limit,” Phys. Rev. B Condens. Matter **49**(24), 17342–17350 (1994). [CrossRef] [PubMed]

**17. **J. P. Huang and K. W. Yu, “Optical nonlinearity enhancement of graded metallic films,” Appl. Phys. Lett. **85**(1), 94–96 (2004). [CrossRef]

**18. **L. F. Zhang, J. P. Huang, and K. W. Yu, “Gradation-controlled electric field distribution in multilayered colloidal crystals,” Appl. Phys. Lett. **92**(9), 091907 (2008). [CrossRef]

**19. **J. Sancho-Parramon, V. Janicki, and H. Zorc, “On the dielectric function tuning of random metal-dielectric nanocomposites for metamaterial applications,” Opt. Express **18**(26), 26915–26928 (2010). [CrossRef] [PubMed]

**20. **M. G. Blaber, M. D. Arnold, and M. J. Ford, “A review of the optical properties of alloys and intermetallics for plasmonics,” J. Phys. Condens. Matter **22**(14), 143201 (2010). [CrossRef] [PubMed]

**21. **See, e.g.,B. Edwards, A. Alu, M. G. Silveirinha, and N. Engheta, “Reflectionless sharp bends and corners in waveguides using epsilon-near-zero effects,” J. Appl. Phys. **105**(4), 044905 (2009) (and references therein). [CrossRef]

**22. **E. O. Liznev, A. V. Dorofeenko, L. Huizhe, A. P. Vinogradov, and S. Zouhdi, “Epsilon-near-zero material as a unique solution to three different approaches to cloaking,” Appl. Phys., A Mater. Sci. Process. **100**(2), 321–325 (2010). [CrossRef]

**23. **K. Maex, M. R. Baklanov, D. Shamiryan, F. Iacopi, S. H. Brongersma, and Z. S. Yanovitskaya, “Low dielectric constant materials for microelectronics,” J. Appl. Phys. **93**(11), 8793–8841 (2003). [CrossRef]

**24. **B. T. Schwartz and R. Piestun, “Total external reflection from metamaterials with ultralow refractive index,” J. Opt. Soc. Am. B **20**(12), 2448–2453 (2003). [CrossRef]

**25. **V. F. Rodríguez-Esquerre, M. Koshiba, H. E. Hernandez-Figueroa, and C. E. Rubio-Mercedes, “Power splitters for waveguides composed by ultralow refractive index metallic nanostructures,” Appl. Phys. Lett. **87**(9), 091101 (2005). [CrossRef]

**26. **T. Tyc and U. Leonhardt, “Transmutation of singularities in optical instruments,” New J. Phys. **10**(11), 115038 (2008). [CrossRef]

**27. **V. G. Kravets, F. Schedin, and A. N. Grigorenko, “Plasmonic blackbody: almost complete absorption of light in nanostructured metallic coatings,” Phys. Rev. B **78**(20), 205405 (2008). [CrossRef]

**28. **J. Hao, J. Wang, X. Liu, W. J. Padilla, L. Zhou, and M. Qiu, “High performance optical absorber based on a plasmonic metamaterial,” Appl. Phys. Lett. **96**(25), 251104 (2010). [CrossRef]

**29. **Z. Fan, R. Kapadia, P. W. Leu, X. Zhang, Y. L. Chueh, K. Takei, K. Yu, A. Jamshidi, A. A. Rathore, D. J. Ruebusch, M. Wu, and A. Javey, “Ordered arrays of dual-diameter nanopillars for maximized optical absorption,” Nano Lett. **10**(10), 3823–3827 (2010). [CrossRef] [PubMed]

**30. **V. G. Kravets, S. Neubeck, A. N. Grigorenko, and A. F. Kravets, “Plasmonic blackbody: strong absorption of light by metal nanoparticles embedded in a dielectric matrix,” Phys. Rev. B **81**(16), 165401 (2010). [CrossRef]

**31. **M. Gaudry, J. Lerme, E. Cottancin, M. Pellarin, J.-L. Vialle, M. Broyer, B. Prevel, M. Treilleux, and P. Melinon, “Optical properties of (Au_{x}Ag_{1-x})_{n} clusters embedded in alumina: Evolution with size and stoichiometry,” Phys. Rev. B **64**(8), 085407 (2001). [CrossRef]

**32. **A. K. Sharma and G. J. Mohr, “On the performance of surface plasmon resonance based fibre optic sensor with different bimetallic nanoparticle alloy combinations,” J. Phys. D Appl. Phys. **41**(5), 055106 (2008). [CrossRef]

**33. **J. C. R. Reis, T. P. Iglesias, G. Douhéret, and M. I. Davis, “The permittivity of thermodynamically ideal liquid mixtures and the excess relative permittivity of binary dielectrics,” Phys. Chem. Chem. Phys. **11**(20), 3977–3986 (2009). [CrossRef] [PubMed]

**34. **C. R. Simovski and S. A. Tretyakov, “Local constitutive parameters of metamaterials from an effective-medium perspective,” Phys. Rev. B **75**(19), 195111 (2007). [CrossRef]

**35. **C. R. Simovski, M. Popov, and S. He, “Dielectric properties of a thin film consisting of a few layers of molecules or particles,” Phys. Rev. B **62**(20), 13718–13730 (2000). [CrossRef]

**36. **C. M. Soukoulis and M. Wegener, “Past achievements and future challenges in the development of three-dimensional photonic metamaterials,” Nat. Photonics **5**, 523–530 (2011).

**37. **J. F. Zhou, T. Koschny, M. Kafesaki, and C. M. Soukoulis, “Negative refractive index response of weakly and strongly coupled optical metamaterials,” Phys. Rev. B **80**(3), 035109 (2009). [CrossRef]

**38. **D. R. Smith, S. Schultz, P. Markoš, and C. M. Soukoulis, “Determination of effective permittivity and permeability of metamaterials from reflection and transmission coefficients,” Phys. Rev. B **65**(19), 195104 (2002). [CrossRef]

**39. **C. H. Gan and P. Lalanne, “Well-confined surface plasmon polaritons for sensing applications in the near-infrared,” Opt. Lett. **35**(4), 610–612 (2010). [CrossRef] [PubMed]

**40. **R. J. Pollard, A. Murphy, W. R. Hendren, P. R. Evans, R. Atkinson, G. A. Wurtz, A. V. Zayats, and V. A. Podolskiy, “Optical nonlocalities and additional waves in epsilon-near-zero metamaterials,” Phys. Rev. Lett. **102**(12), 127405 (2009). [CrossRef] [PubMed]

**41. **M. G. Blaber, M. D. Arnold, and M. J. Ford, “Designing materials for plasmonic systems: the alkali-noble intermetallics,” J. Phys. Condens. Matter **22**(9), 095501 (2010). [CrossRef] [PubMed]

**42. **S. K. Golden and G. Papanicolaou, “Bounds for effective parameters of heterogeneous media by analytic continuation,” Commun. Math. Phys. **90**(4), 473–491 (1983) (and references therein). [CrossRef]

**43. **E. Tuncer, “Extracting the spectral density function of a binary composite without a priori assumptions,” Phys. Rev. B **71**(1), 012101 (2005). [CrossRef]

**44. **E. Tuncer, “Geometrical description in binary composites and spectral density representation,” Materials **3**(1), 585–613 (2010). [CrossRef]

**45. **D. Zhang and E. Cherkaev, “Pade approximations for identification of air bubble volume from temperature- or frequency-dependent permittivity of a two-component mixture,” Inv. Probl. Sci. Eng. **16**(4), 425–445 (2008). [CrossRef]

**46. **C. Bonifasi-Lista and E. Cherkaev, “Electrical impedance spectroscopy as a potential tool for recovering bone porosity,” Phys. Med. Biol. **54**(10), 3063–3082 (2009). [CrossRef] [PubMed]

**47. **A. V. Goncharenko, V. Lozovski, and E. F. Venger, “Effective dielectric response of a shape-distributed particle system,” J. Phys. Condens. Matter **13**(35), 8217–8234 (2001). [CrossRef]

**48. **E. Cherkaev and M.-J. Y. Ou, “Dehomogenization: reconstruction of moments of the spectral measure of the composite,” Inv. Probl. **24**(6), 065008 (2008). [CrossRef]

**49. **G. A. Niklasson and C. G. Granqvist, “Optical properties and solar selectivity of coevaporated Co-Al_{2}O_{3} composite films,” J. Appl. Phys. **55**(9), 3382–3410 (1984). [CrossRef]

**50. **S. Riikonen, I. Romero, and F. J. Garcia de Abajo, “Plasmon tunability in metallodielectric metamaterials,” Phys. Rev. B **71**(23), 235104 (2005). [CrossRef]

**51. **P. Mallet, C. A. Guérin, and A. Sentenac, “Maxwell-Garnett mixing rule in the presence of multiple scattering: Derivation and accuracy,” Phys. Rev. B **72**(1), 014205 (2005). [CrossRef]

**52. **C. J. F. Böttger, *Theory of Electric Polarization* (Elsevier, 1952).

**53. **A. N. Lagarkov and V. N. Kisel, “Losses in metamaterials: Restrictions and benefits,” Physica B **405**(14), 2925–2929 (2010). [CrossRef]

**54. **S. Xiao, V. P. Drachev, A. V. Kildishev, X. Ni, U. K. Chettiar, H. K. Yuan, and V. M. Shalaev, “Loss-free and active optical negative-index metamaterials,” Nature **466**(7307), 735–738 (2010). [CrossRef] [PubMed]

**55. **S. M. Anlage, “The physics and applications of superconducting metamaterials,” J. Opt. **13**(2), 024001 (2011). [CrossRef]