1. Introduction

Homogenization of inhomogeneous materials (i.e. finding their effective conductivity, permittivity, permeability, elastic moduli, etc.) has long been an important problem in science and engineering. Recent interest in metamaterial (MM) design has given renewed impetus to this activity. Despite a widespread interest in this topic, currently only a few closed-form analytical results are known which may be applied to MMs of special kind [1]. Among them, the most known are those describing one-dimensional (1d) systems, originally introduced by Voigt and Reuss for effective elastic moduli [2, 3]. For two-dimensional (2d) MMs, there are also few microgeometries which yield explicit formulae for their effective moduli (see, e.g., [1, 4]). A rather general and accurate approach for homogenization of 2d periodic media has been proposed by Krokhin et al. [5]. However, it inlolves dealing with matrices of big size, requires a lot of computing resources and, in addition, can suffer from slow convergence.

As was noticed by Sahimi [6], systems which are usually referred to as 1d, are in fact quasi-1d systems. A logical question arises: what if there is a small deviation from one-dimensionality, i.e., microgeometry is almost 1d? From an intuitive viewpoint, an analytical solution should exist in this case, too. One way to address this issue is to invoke the formalism of the spectral density function [7]. In this paper, we show how to find an approximate solution to the problem involving phase content distribution function. For definiteness, in what follows we consider the problem in terms of the effective permittivity, but other effective moduli can be also considered. Besides, here we deal with homogenization in the quasistatic limit. Such an approximation is thought to be valid when a scale size for spatial variation of the structure under study is much less than the wavelength. The effects of spatial dispersion are beyond the scope of our work.

In fact, MMs with 1d-microgeometry are of interest by themselves offering a lot of unusual phenomena and potential applications, such as, for example, etraordinary transmission [8], optical filters, mirrors, and polarizers [9, 10], biosensing [11], polarization conversion [12], or subwavelength imaging [13–15]. Of special interest among them are the so-called hyperbolic MMs, known also as indefinite media (see, e.g., [14, 16–18]). In turn, quasi-1d MMs, introducing additional degrees of freedom, could offer even more possibilities. So, they could help in the design of structures possessing filtering properties, ultralow refractive index (ULRI) and enhanced absorption over a wavelength range [19] or having broadband low permittivity (broadband epsilon-near-zero (ENZ) MMs) [20, 21]. In our previous studies [19–21], we considered quasi-1d MMs with a complex unit cell and determined microgeometrical parameters of the cell solving an inverse problem. However, if the number of the parameters is large, solving the inverse problem can present difficulties. Besides, there we dealt only with one component of the effective permittivity tensor. In the present work, we follow a different route. It consists in setting simple functional forms for the cell microgeometry and finding both diagonal components of the tensor. In the considered examples, the effective permittivity may be written in closed form that simplifies its analysis from the point of view of potential applications. So, we show that ENZ, ULRI, and high refractive index (HRI) MMs with extended bandwidth can be designed within the framework of our approach.

2. Basic formalism

Let us first consider the case of an ideal one-dimensional MM. The simplest 1d geometry consists of alternating layers with different permittivities. For definiteness, we deal with a material, which is periodic in the y-direction, and homogeneous in the xz-plane. So, the y-axis is the direction along which the local permittivity can change. To validate the efeective permittivity formalism, the period S_y is considered to be much less than the wavelength (quasistatic approximation). Then, the effective permittivity can be represented as a tensor of the form ε̃ = diag(ε_x, ε_y) with

In the quasistatic approximation, Eqs. (1) and (2) are exact; the only limitation is that they do not take into account spatial dispersion. Let us now allow a deviation from one-dimensionality. To do this, we turn to specific microgeometry represented in Fig. 1. Here we show the unit cell of a 2d MM made of two materials with the permittivities ε₁ and ε₂. In the following, the phases 1 and 2 are called the slab and the host, respectively. Therefore, the above MM is a periodic array of parallel slabs of varied width. Then, the local width of the slab normalized to the size of the unit cell along the x–axis S_x can be defined as the local filling factor of the phase 1, f(y). Obviously, the total filling factor of the phase 1 is $f_1={\int }_0^{S_y}f(y)dy$ and that of the phase 2 is f₂ = 1 − f₁. Next we break the unit cell down into thin layers of the thickness dy (one such a layer is shown in Fig. 1). The assumption, used in the following analysis, is that the interface between two phases remains almost parallel to the y–axis within each layer. To meet this condition, the function f(y) should, obviously, be smooth and satisfy the inequality

 

Fig. 1 Schematic of the unit cell representing geometry under consideration. The arrows on the right side of the figure show possible orientations of the electric field.

Download Full Size | PDF

The accuracy of Eqs. (6) and (7) is controlled by In eq. (3). Numerical estimations of the accuracy of Eq. (7) for rectangular unit cells with different aspect ratios are given in [20, 21]. Equation (7) was earlier obtained in [21] and, in a slightly different form, in [19]. There, imposing some constraints on ε_y, we solved an inverse problem and found the piece-wise function f. Such an approach can, however, involve computational difficulties. Here, we suggest to parametrize the distribution function f(y) using only a few parameters that considerably simplifies the practical use of Eqs. (6) and (7).

3. Particular cases of the distribution function f(y′)

While there is a lot of options for the choice of the distribution function, in the following we restrict ourselves to the two simplest options. Those are the functions in the form of the linear and parabolic profile (see Fig. 2). Both functions have three free parameters: y′₀, f_min, and f_max. The linear profile can be analytically represented as

Substituting Eq. (9) into Eqs. (6) and (7), one obtains

Analogously, substituting Eq. (10) into Eqs. (6) and (7), one has

 

Fig. 2 The linear and parabolic profiles of the metal phase distribution along the y-axis.

Download Full Size | PDF

4. Examples of application: Numerical results

The operational bandwidth is a critical issue for many optical devices. At the same time, most of reported MMs retain their useful properties only within a relatively narrow frequency band, and extending the bandwidth seems to be a challenging task. A possible solution to this problem with the use of the quasi-1d MMs has now been demonstrated for ENZ, ULRI, and HRI media. In all examples given below we take Teflon with the refractive index $n_2\equiv {\epsilon }_2^{1/2}=1.35$ as the dielectric phase and gold as the metal phase. The permittivity of gold has been taken with the use of the Drude-critical points model as suggested in [23].

4.1. ENZ metamaterials

Earlier, it has been shown that just the ε_y component of the dielectric tensor can be easily nulled or kept small over a broad frequency band [20, 21, 24, 25]. This is related to the feature of multilayered MMs that can be formulated as follows: if the electric field is oriented normally to the layers and at least one of the layer permittivities is zero, the corresponding component of the effective permittivity tensor must be zero, too. Then, having the parallel layers with the permittivities nulled at closely located frequencies, it becomes possible to bring the effective permittivity curve closer to zero within an entire frequency range. Although the above statement is correct (in the quasistatic approximation) for ideal 1d MMs only, much of the same is true of quasi-1d MMs. Thus, mathematically the problem can be reduced to minimizing the norm ||Reε_y[ε₁(ω), ε₂, f(y′₀, f_min, f_max)]|| on a frequency band [ω₁, ω₂] with obvious constraints 0 ≤ y′₀ ≤ 1/2 and 0 ≤ f_min ≤ f_max ≤ 1, that allows one to find the needed geometrical parameters.

It is necessary to keep in mind that a minimization procedure may yield several solutions. As we have made sure, for the chosen materials and frequency range (the visible and near IR), at least two solutions occur. One of them is shown in Fig. 3, where we present the results of Reε_y fitting to zero for the bands 570 – 670 nm, 570 – 720 nm, and 570 – 770 nm (here we use the parabolic profile for the distribution function; the use of the linear profile yields worse fitting). The root mean square (rms) of the deviation of the fit from zero for each band is 8.8 · 10⁻⁴, 0.0015, and 0.0022, respectively. It should be also noted that the shift of the actual band to the long wavelength side allows one to considerably improve the accuracy of fitting. For example, rms = 1.3 · 10⁻⁴ for the band 700 – 800 nm (see Fig. 4, blue solid curve). In this figure, two solutions for Reε_y, fitted to zero over the band 700 – 800 nm, as well as the corresponding imaginary parts of the effective permittivity are shown.

 

Fig. 3 The real part of the effective permittivity of designed MMs fitted to zero over the bands of 570 – 670 nm (1), 570 – 720 nm (2), and 570 – 770 nm (3). The fitted parameters are: y′₀ = 0.248, f_min = 0.043, f_max = 0.627(1); y′₀ = 0.25, f_min = 0.035, f_max = 0.668(2);y′₀ = 0.252, f_min = 0.03, f_max = 0.699(3).

Download Full Size | PDF

 

Fig. 4 The effective permittivity of two designed MMs with Reε_y fitted to zero over the band of 700 – 800 nm. The fitted parameters are: y′₀ = 0.259, f_min = 0.021, f_max = 0.975 (blue curve) and y′₀ = 0.177, f_min = 0.063, f_max = 0.103 (red curve).

Download Full Size | PDF

4.2. ULRI metamaterials

The problem of the design of broadband MMs with the ultra-low refractive index n^* is similar to that for ENZ MMs [19]. In this case, it can be reduced to minimizing the norm ||Re{ε_y[ε₁(ω), ε₂, f(y′₀, f_min, f_max)]}^1/2 − n^*|| on a frequency band [ω₁, ω₂] [19, 21]. When minimizing the norm and finding the sought parameters y′₀, f_min, and f_max, multiple solution branches can coexist for the same n^*. As an example, in Fig. 5 we show the effective refractive index of the designed MM fitted to n^* = 0.25 over the band λ = 620 – 720 nm (the red edge of the visible spectrum). As before, we have used the parabolic profile for the distribution function, although the use of the linear profile also yields satisfactory results. The four curves (1–4) correspond to four sets of the parameters y′₀, f_min, and f_max which allow one to obtain local minima of the above norm. The root mean square of the deviation of the fit from n^* for each solution branch is shown in Fig. 6. As one can see, only one solution branch (curve 2) exists for n^* < 0.1475, while another branch (curve 4) provides appropriate solution (with small rms) as n^* → 1. At n^* = 0.16, the minimal rms can be as small as 0.0045. At the same time, the root mean square rapidly rises as n^* becomes smaller than 0.14.

 

Fig. 5 The effective refractive index of designed MMs fitted to n^* = 0.25 over the band of 620 – 720 nm. The fitted parameters are: y′₀ = 0.492, f_min = 0, f_max = 0.769 (1); y′₀ = 0.498, f_min = 0.9425, f_max = 1 (2); y′₀ = 0.396, f_min = 0, f_max = 0.1615 (3); y′₀ = 0.436, f_min = 0.09, f_max = 0.15 (4).

Download Full Size | PDF

 

Fig. 6 The four branches for the root mean square curves vs n^*.

Download Full Size | PDF

4.3. HRI metamaterials

The straightforward way to design broadband MMs with high refractive index is to minimize the norm ||Re{ε_x[ε₁(ω), ε₂, f(y′₀, f_min, f_max)]}^1/2 − n^*|| on a frequency band [ω₁, ω₂]. The effective refractive index for two possible solutions (both are for the parabolic distribution function) for n^* = 4.25 are shown in Fig. 7. We note that although both curves have been obtained for the same actual band 620 – 720 nm, one of the solutions (solid curve) provides appropriate fitting over a broader band. Furthermore, as can be seen from Fig. 8, three solution branches can coexist below n^* = 3 and two branches above that. One of them (curve 1) provides lower rms below n^* = 4.25, while another provides lower rms above it. For example, at n^* = 4.5 the best fit provides rms ≃ 0.02 and then the root mean square rapidly rises with n^*. It should be also noted that shifting the actual band to the long wavelength side allows one to considerably extend the bandwidth providing high values of the effective refractive index and keeping rms relatively small. For example, fitting Re(ε_x)^1/2 to n^* = 5 over the band of 660 – 860 nm yields rms = 0.006, while fitting to n^* = 5.5 over the band of 700 – 900 nm yields rms = 0.0091 (see Fig. 9). Fitting to n^* = 6 over the band of 770 – 1020 nm (not shown here) yields the root mean square as small as 0.0116.

 

Fig. 7 The effective refractive index of designed MMs fitted to n^* = 4.25 over the band of 620 – 720 nm. The fitted parameters are: y′₀ = 0.063, f_min = 0.817, f_max = 0.936 (solid line); y′₀ = 0.065, f_min = 0.805, f_max = 0.918 (dashed line). The root mean square is about 0.017 for both curves.

Download Full Size | PDF

 

Fig. 8 The three branches for the root mean square curves vs n^* in the HRI regime.

Download Full Size | PDF

 

Fig. 9 The effective refractive index of designed MMs fitted to n^* = 5 over the band of 660 – 860 nm (rms = 0.006) and to n^* = 5.5 over the band of 700 – 900 nm (rms = 0.0091). The fitted parameters are: y′₀ = 0.042, f_min = 0.856, f_max = 0.956 (solid line); y′₀ = 0.041, f_min = 0.87, f_max = 0.953 (dashed line).

Download Full Size | PDF

5. Discussion

Why quasi-1d MMs are promising for the design of ULRI ans HRI media with extended bandwidth? To address this point, it seems useful to consider what happens in the narrowband case. For any two-phase composite, the diagonal components of the effective permittivity tensor must lie within a region in the complex plane (Reε, Imε) which is bounded by the straight line, connecting ε₁ and ε₂, ε_u = f₁ε₁ + f₂ε₂, and by the arc, connecting the origin, ε₁, and ε₂, ε_l = (f₁/ε₁ + f₂/ε₂)⁻¹ (the above bounds are known as upper and lower Wiener’s bound, respectively, in terms of the effective conductivity). In turn, the allowable values of the effective refractive index must lie within a region in the complex plane (Ren, Imn), bounded by the curves $n_u={\epsilon }_u^{1/2}$ and $n_l={\epsilon }_l^{1/2}$. The formal similarity between ENZ and ULRI metamaterials consists in the fact that both (ENZ and ULRI) regimes are realized close to upper Wiener’s bound. Above we have noticed why this takes place for the ENZ regime. As to the ULRI regime, this takes place because, as one can prove, always Ren_u ≤Ren_l, and no allowed value of ε_eff exists such that $\text{Re}{\epsilon }_{\mathit{eff}}^{1/2}<\text{Re}{n}_u$.

The function

Looking at Eqs. (15),(16), and (18), it is easy to see that reducing n_min can be achieved by reducing ε₂. Besides, if the metal permittivity is described by Drude’s model, at low frequencies n_min is large: $n_{\mathit{min}}\approx \frac{\gamma }{\omega }\sqrt{{\epsilon }_2}\to \infty$ as ω → 0 (here γ is the damping coefficient or collision frequency). On the other hand, it is also large at high frequencies, because the dielectric constrast Δ is relatively small in this case. Thus, at some frequency, n_min takes its minimum value.

Because, as was noticed above, Ren_u ≤Ren_l, the point with the maximum value of the effective refractive index lies on the lower Wiener’s bound. That is why just Eq. (6) (which, as we remind, corresponds to electric field oriented normally to the wires) is convenient to realize the HRI regime. Upon differentiating Ren_l(f) and equating the result to zero, it is possible to find $f_{\mathit{max}}^{*}$ at which the real part of the effective refractive index takes its maximum value. After simple but bulky calculations, the final result may be written as follows:

The dependencies n_min(λ) and n_max(λ) are shown in Fig. 10. So, n_min is large at high frequencies, reduces with reducing the frequency, and takes its smallest value (0.078) at about λ = 860 nm. Then, it slowly rises with λ, i.e., with a decrease in the frequency. The results of calculations performed with the use of exact and approximate formulas are almost indistinguishable at wavelengths above λ = 560 nm. At the same time, n_max monotonically rises with λ.

 

Fig. 10 The frequency dependencies of $n_{\mathit{min}}={\mathit{Ren}}_u(f_{\mathit{min}}^{*})$ calculated with the use of the exact Eqs. (15) and (16) (solid curve), with the use of the approximate Eq. (18) (dashed curve), and $n_{\mathit{max}}={\mathit{Ren}}_l(f_{\mathit{max}}^{*})$ (dotted curve).

Download Full Size | PDF

The above results for n_min and n_max have been obtained for a fixed frequency (the narrowband case). Obviously, n_min cannot be reduced, while n_max cannot be enlarged when extending the bandwidth.

The presence of several solution branches poses a question which one is more preferable in practice. At first glance, it makes sense to choose a solution for which rms is minimal. However, this is not always so. Let us rewrite the condition of quasi-one-dimensionality, In eq. (3), in the form

In the above consideration, we have tailored the real parts of both the effective permittivity and refractive index; nothing has been said about their imaginary parts being resulted from losses. Usually, high losses prevent the practical realization of the unique applications of MMs; sometimes, however, they are desirable [26, 27]. As can be seen from Fig. 4, the imaginary part of the effective permittivity can differ markedly for designed ENZ MMs (in particular, this depends on the total content of dissipative (metal) phase). So, it can take rather large values (dashed blue curve), as well as moderate values (dashed red curve). Behavior of the imaginary part of the effective refractive index can be inferred from Fig. 11 where we show both lower and upper Wiener bounds in the complex (Ren,Imn) plane for three wavelengths, 620, 670, and 720 nm. We note that in the ULRI regime, which is realized near the upper Wiener bound, the imaginary part of the effective refractive index distinctly rises as n^* → n_min. At the same time, it can be small at moderate values on n^* and very small when n^* is about unity or slightly above it (this is because the ULRI, as well as ENZ regime are nonresonant in character [19, 20]). As to the HRI regime, which is realized near the lower Wiener bound, the imaginary part of the effective refractive index is always high (obviously, this takes place since this regime, in contrast, is resonant in character).

 

Fig. 11 The upper (solid curves) and lower (dashed curves) Wiener’s bounds in the complex (Ren, Imn) plane for λ = 620, 670, and 720 nm.

Download Full Size | PDF

6. Generalizations and limitations

In the above consideration, we dealt with parallel slabs. At the same time, the 3d generalization of the problem is straightforward. Indeed, in this case the effective permittivity can be considered as the tensor ε̃ = diag(ε_x, ε_z, ε_y) where the ε_y component is of the same form as in Eq. (2),

How accurate are Eqs. (6) and (7)? This question is related to the applicability of the used homogenization technique; of principal importance in this connection are the effects of nonlocality (spatial dispersion) neglected here.

One can distinguish two lines of research on nonlocal effects in MMs. The first line takes into account intrinsic nonlocality of constituents (material nonlocal response). For metals, nonlocal effects result from the quantum nature of the free electrons and their interactions (longitudinal pressure waves of the free-electron plasma [29]). For noble metals, it is generally agreed that such nonlocal effects become significant when characteristic device dimensions approach 1 – 10 nm. As has been recently shown [30], in the ENZ regime there are no evidences that nonlocal response can cast doubt on the regular local-response approximation. In addition, perfect imaging has only a weak dependence on intrinsic metal nonlocality.

The second line of research deals with materials’ response described within the framework of the local-response approximation. In this case, the need to take into account the spatial dispersion and associated magnetoelectric coupling [31] for the accurate description of the effective parameters of MMs has been noted in many papers (see, e.g., [32] and references therein). One of the reasons for strong spatial dispersion is resulted from excitation and coupling of surface plasmon polaritons at the metal-dielectric interfaces [33, 34]. In some cases, e.g., for perfectly conducting parallel thin wires, spatial dispersion is strong and cannot be neglected [35]. In fact, this problem in not trivial and its detailed discussion does not enter into the scope of the present paper. It should be, however, noted that, to a first approximation, it looks reasonable to consider the spatial dispersion effects inherent in quasi-1d MMs in terms of those inherent in ideal 1d MMs. At the same time, for (ideal) 1d MMs the similar problem has been considered in a number of papers (see [15] and references therein). In particular, it has been shown that the spatial dispersion can be suppressed in the near-infrared and visible ranges, where the kinetic inductance of metal nanowires can strongly exceed their magneto-static inductance. The vanishingly small spatial dispersion in the visible and near-infrared ranges for an assembly of parallel gold nanorods grown into a porous alumina template has been demostrated in experiment including both ENZ and ULRI regimes [11]. Another recipe to suppress the spatial dispersion in wired media is to load the wire medium with metal patches in order to increase the effective capacitance of the wires [36]. Furthermore, because the waves propagating across the layers cannot couple to surface plasmon polaritons, the spatial dispersion can occur only for oblique propagation and should be weak for on-plane propagation (at k_y = 0); in addition, its strength can be controlled by changing a ratio between the thicknesses of metal and dielectric layers [33].

Some data evidence that even in the near-infrared range the magnetic permeability of metal wire arrays can differ from unity close to the magnetic antiresonance region [37]. Obviously, in such cases our results for the ULRI and HRI MMs should be reconsidered. When taking into account the magnetic properties, we note also that new interesting phenomena like resonant absorption and field enhancement can occur in the so-called transition MMs near the zero-index point or near two spatially separated zeros of the permittivity and permeability [38].

In our calculations, we used bulk values of the permittivity of gold. Meantime, for thin metal layers it can differ from the bulk permittivity. So, the damping constant γ should be corrected to take into account the interfacial scattering effects, and may be written as γ = γ_b + Av_F/t where γ_b is the bulk damping constant, A is a constant which depends on the surrface scattering mechanism, v_F is the Fermi velocity of electrons, and t is the characteristic thickness of the layer. However, as we noticed earlier [20], the effect of the surface scattering can be neglected when the electric field is oriented along the interfaces, i.e., for the ε_y component of the dielectric tensor. At the same time, this effect can take place as the electric field is perpendicular to the layers which, in turn, are thin enough. In this case the metal permittivity ε₁ can change, that has to be taken into account in practical calculations.

7. Conclusion

To conclude, in this paper we have considered two principal issues. First, we have pointed the way to homogenization of MMs made of parallel slabs, when MM geometry is not too different from the 1d geometry. Two analytical expressions for two diagonal components of the effective permittivity tensor have been written in terms of the coordinate-dependent distribution function of the phase content, and then two simpliest examples of this function, which set the linear and parabolic profiles of the slab thickness, have been considered. Second, we have demonstrated that such MMs hold promise to extend the bandwidth of media with very small (close to zero) effective permittivity, as well as ultra-low (less than unity) and extraordinary high effective refractive index. In addition, generalization of our method for the 3d case has been outlined, and limitations, such as the spatial dispersion and interfacial scattering effects, have been briefly examined.