# Extension of the Multipole Approach to Random Metamaterials

Authors: A. Chipouline; S. Sugavanam; J. Petschulat; T. Pertsch
Journal: Advances in OptoElectronics
 (2012)

---

## Abstract

Influence of the short-range lateral disorder in the meta-atoms positioning on the effective parameters of the metamaterials is investigated theoretically using the multipole approach. Random variation of the near field quasi-static interaction between metaatoms in form of double wires is shown to be the reason for the effective permittivity and permeability changes. The obtained analytical results are compared with the known experimental ones.


---

## Body

## 1. Introduction

Metamaterials are artificial media that allow tailoring the macroscopic properties of the light propagation by a careful choice of a design for the microscopic unit cell (called the meta-atom). By controlling the geometrical shape and the material dispersion of the meta-atom, novel effects such as negative refraction [1–3], optical cloaking [4–9], as well as series of optical analogues to phenomena known from different disciplines in physics could be observed [10–14]. Despite the possibility to rely on rigorous computations for describing the light propagation on the microscopic level, an enduring problem in metamaterial research is the question on how the effective material tensor looks like for a certain metamaterial. A simple and versatile analytical model describing propagation of electromagnetic waves in metamaterials has been recently developed following classical approach of Maxwell equation averaging procedure [15]. This transition from the microscopic to macroscopic system of Maxwell equations takes into account all peculiarities of carriers dynamics under the action of the resulted electromagnetic field through the introduction of multipole moments which are supposed to be represented as the functions of the macroscopic electric and magnetic fields [16]. One of the great advantages of this model is the ability to evaluate straightforwardly influence of the charge dynamics of the meta-atoms on the effective properties of the metamaterials. In fact, the multipole moments are calculated through the averaged charge dynamics in the meta-atoms. Any factors influencing the charge dynamics (e.g., interaction between the meta-atoms, extra coupling of the meta-atoms with the other objects, etc.) cause the changes in the multipole expressions, which in turn change the effective parameters. It is important for the analysis presented here that the interaction between the meta-atoms and hence its influence on the effective permittivity and permeability can be straightforwardly taken into account [17].The interaction between the small particles, both dielectric and metallic, and propagation of an optical excitation in a regular chain of such particles have been intensively investigated [18–23]. Interest to the chains of metallic nanoparticles is stipulated mainly by the request for subwavelength guiding structures for a new generation of the optoelectronic components for communication and information processing. Nevertheless, theoretical tools for the modeling of these chains (irrespective to the nature and sizes) remain invariant: the electromagnetic excitation in the particles are supposed to be described by taking into account all possible eigen modes [18, 20] and interaction between all particles in a chain. There are several approximations which are typically accepted in this kind of problems. First, depending on the size of the particles, the model can be restricted by consideration of dipole moment only (for metallic nanoparticles) [19, 23]; the higher moments can be taken into consideration as well in the case of investigation of magnetic response [22, 24]. Usually, for the problem of electromagnetic excitation propagation along the chain, the dipole approximation is enough [25], provided distance between particles is not less than about three times their sizes. Second, the interaction between the particles in a chain can be considered in the frame of the quasi-static approximation, where no retardation between particles is retained; otherwise, interaction between dipoles contains terms proportional to the 1/r and 1/r2 in addition to the quasi-static term 1/r3 (r is the distance between dipoles). The problem possesses an exact solution for the infinite chain in the quasistatic case, while taking into consideration the retardation leads to known math difficulties and requires continuation into the lower half frequency plane [19]. Consideration of the finite chain is free from these excessive math problems but can be treated only numerically; the respective solutions for both longitudinal and transverse modes are presented in [19, 26].Natural expansion of the developed models of the electromagnetic excitation transport on a chain of particles with randomly varying parameters revealed several interesting peculiarities. The problem of wave propagation through disordered systems attracts great attention in both quantum and classical physics [27]. In disordered chains of different dimensions, destructive interference between scattered waves gives rise to an existence of the localized modes, exponentially decaying in space—this effect has been originally found in solid state physics and is known as Anderson localization [28]. The existence of delocalized modes that can extend over the sample via multiple resonances and have a transmission close to 1 was found in [29, 30] and experimentally confirmed in [31, 32]. Disorder-induced change of the guiding properties in a chain of plasmonic nanoparticles under small random uncontrollable disorder was considered in [26], and analogy of the Anderson localization in a chain of such particles was theoretically investigated in [33]. In the present analysis, the effect of Anderson localization is not considered; nevertheless, it is believed, that the developed analytical tool in this paper turns out to be suitable for the treatment of the similar effect in metamaterials with different types of disorder.An influence of various types of disorder on the effective properties of the metamaterials was intensively investigated as well. Light propagation and Anderson localization in superlattices were theoretically considered in [34, 35] using the model of multilayered system with phenomenological permittivity and permeability (positive and negative) in each of the layers. The effect of the statistical distribution of the sizes of the meta-atoms on the increase of losses in the operation frequency band was considered in [36] using generalized Clausius-Mossotti relation. A significant influence of a small (10%) deviation of the parameters of the microscopic resonances on the propagation wave in a wide frequency range was found in [37] using quasistatic expressions for the effective parameters. Averaging of the Lorenz-type expressions for the effective permittivity and permeability using a phenomenological probability distribution function showed that passband and negative refraction are still present under small positional disorder [38]; the results were proven experimentally as well. Interaction in a chain of magnetic particles and its influence on the effective permeability were investigated in [24]. Using the introduced concept of “coherent” and “incoherent” metamaterials, authors of [39] showed that the influence of disorder on long-range correlated metamaterials is significantly more pronounced in comparison with the same effect in short-range ordered metamaterials. Random variation of interaction between meta-atoms was shown to be a main reason for the disappearance of the long-range correlation and consequently of the “coherent” state [39].In the presented work, attention is primarily devoted to the extension of the multipole approach to describe in-plane disorder in metamaterials, which means the randomness in positioning of the meta-atoms in the plane of the substrate; see Figure1. In [26], it was shown that the variations of the electromagnetic properties of the inclusions are less important than the disorder in their positions. Metamaterials formed by a self-organization display exactly this kind of disorder [40–43]. Results of control experiments with the 2D metamaterials exhibiting such in-plane disorder [44] are used as a test of the model. The most notable discovery is the fact that although disorder has a deterrent effect on the permittivity, the permeability seems to remain practically unaffected. A theoretical model for such a class of random metamaterials should reproduce these observations.(a) Regular and (b) laterally (alongz direction) random positioning of the meta-atoms in a metamaterial. Plane wave propagation is in y direction, and electric field polarization coincides with the x axis and with the elongation of the nano wires. Note that only one layer, which the metamaterial consists of, is shown in the figure.
(a)
(b)The qualitative explanation of the influence of the spatial disorder on the effective parameters is in following. The positional disorder creates different conditions for the charged dynamics in the meta-atoms due to the interaction between them [21, 22]. This in turn leads to the changes of the averaged dipole, quadrupole, and magnetic dipole moments of the media and results in the changes of the effective parameters, which are expressed through these averaged multipole moments [15, 45].This qualitative hypothesis requires further development of the existing theoretical multipole model; in particular, the interaction between the meta-atoms [17] has to be incorporated and adopted to the random character of this interaction. Let us assume that the charge dynamics in the microscopic multipole moments of the meta-atoms depends on the distance δk between them (see Figure 2). Following the approach of [15], it is necessary to average the resulted charge dynamics in the multipole moments over all possible representations; in other words, the microscopic multipole moments have to be additionally averaged over all possible distances between the meta-atoms, which mathematically is expressed as an integral over a probability distribution function PDF(δk); namely,
(1)χmacro=∫χmicro(δk)PDF(δk)dδk,
Here, χmicro(δ) is the microscopic multipole moment of the meta-atoms, and PDF(δk) governs the distribution over all possible interseparation distances δk in a randomly arranged ensemble of the meta-atoms. In case of regular spatial distribution, each metaatom is affected by the same fields, PDF(δk) is reduced to the delta function, and averaging (1) restores the microscopic multipole moments.Figure 2
Geometry for the probability function elaboration, the spheres show meta-atoms. The first row shows a regular arrangement of the meta-atoms, where each meta-atom occupies the center of a slot of length equal to the mean period. The second row depicts an arrangement of the meta-atoms exhibiting random uncorrelated positional disorder (denoted byρk), the extent of the disorder being governed by the PDF(ρk) as shown in the last row. The interseparation δk between the two subsequent meta-atoms is a function of the random variables ρk and ρk-1, and the analytic form of the PDF(δk) can be obtained by the use of the statistical methods if the analytic form of the PDF(ρk)  is given.The quest to obtain suchPDF(ρk) and the effort to incorporate the effect of disorder into the existing multipole model are discussed in details here. The paper is organized as follows. First, the probability model used to incorporate positional disorder into the multipole theory is described. As a test of principle and in order to create a systematic model, the approach is then applied to the simple case of randomly arranged dipoles. Then the treatment is extended to the case of randomly arranged quadrupoles. The probabilistic approach is applied to the specific case of randomly positioned meta-atoms, and the obtained results are compared with the experimental observations [44]; the mathematical procedures used to account for the other forms of disorder are highlighted.The main difference of the approach here in comparison with the previous ones in the use of the multipole model is that the charge dynamics in meta-atoms is primarily considered and calculated taking into account the interaction between meta-atoms, which is expressed as a function of distance between them. Finally, averaging over all possible realization of the inter-metaatom distance gives the expression for the effective parameters. This paper is primarily devoted to the elaboration of the model and to the effective parameters calculation; further applications of the presented approach (disorder in propagation direction, transition “coherent-incoherent” states, influence of the Anderson localization on the effective parameters, etc.) will be done elsewhere. Interaction between meta-atoms is taken into account using a simplest way of dipole-dipole near field interaction in quasi-static limit; extrapolation of the model on the dynamic case is left for the future work. The interaction between quadrupoles is treated the same way, which makes the presented approach suitable for consideration of the magnetic properties of the metamaterials. In spite of the excessive simplification of the interaction, the presented model treats the effective parameters (especially magnetic response) in much more correct way than it was done before by just introduction of permeability and/or magnetic susceptibility, and is believed to provide a suitable platform for analytical or semianalytical treatment of the problems, appearing in the case of disordered metamaterials.


## 2. Modeling of Positional Disorder

The problem of the positional disorder modeling can be tackled in several ways. The most general formulation of the problem requires a Markovian treatment. As an illustration, the one-dimensional equivalent of the problem is considered here. Supposing that the meta-atoms are introduced one by one on a line of given length, the probability that a particle will take up a certain position on the line, and hence the probability of a particular inter separation distance, depends not only on the last particle, but also on the history and existing configuration. This is the essence of the Markovian approach. Standard techniques exist for tackling such problems, formulating a rigorous treatment. Nevertheless, due to its complex nature, an alternative simpler math treatment is developed in the present paper.The math treatment accepted in the presented work is stipulated by the real technological chains for the producing of the nanostructure. When performing control experiments, masks for random metamaterials are manufactured by e-beam lithography methods. The writing algorithms are modified so that instead of a periodic grid a randomized one is generated by the scanner. The extent of randomness can be specifically controlled, and statistically relevant parameters such as the mean period and the variance can be assigned to each mask. Translating the above approach to mathematics, one assumes that the meta-atoms are initially placed in well-defined slots of equal length (see Figure1) and then are perturbed from their mean positions. The perturbation is described by  PDF(ρk). This PDF describes the extent of the perturbation and also ensures that the displacement is not beyond a certain space slot. However, PDF(ρk) is the positional disorder distribution function not the PDF for the metaatom interseparation PDF(δk); the latter has to be found based on given PDF(ρk), and  PDF(δk) can then be used in averaging procedure (1) to obtain the required material parameters. In this work, interseparation PDF(δk)  is obtained by employing a characteristic function approach. The characteristic function Qρ(ω)   is given by a Fourier transformation of given PDF(ρk):
(2)Qρ,k(ω)=∫-∞∞PDF(ρk)eiωtdρk,
where PDF(ρk) satisfies the normalization condition:
(3)∫-∞∞PDF(ρk)dρk=1.
Alternatively, the characteristic function can be considered as an expectation value of function eiωt:
(4)Qρ,k(ω)=〈eiωt〉|PDF(ρk).
The one-dimensional equivalent of the problem is formulated as follows. A periodic arrangement of N meta-atoms on a given length L is considered (see Figure 1). The spacing period for the slots is given by z0=L/N, and the location of the kth metaatom can thus be given as zk=(z0/2)(2k+1). Now, the perturbation of the kth metaatom from its mean position can be given using a random function ρk, such that.
(5)zk=z02(2k+1)+ρk,-zk2<ρk<zk2.
So, the spacing between the meta-atoms is given by.
(6)Δk=zk+1-zk=z0+(ρk+1-ρk)=z0+δk.
The random functions ρk and ρk+1 are completely independent from each other. For the present problem, the random function of interest is
(7)δk=ρk+1-ρk.
The mathematical form of  PDF(δk) has to be found. The form of the above probability function can be obtained by using characteristic functions. The characteristic function formed by a sum or difference of two or more PDFs is nothing but the product of the characteristic functions of the ingredient PDFs. That is, if one is interested in a probability distribution function for variable z, given by
(8)z=y1+y2+⋯+yn,
where y1, y2 and so forth  are independent random functions and  Qi(ω) are their mutually independent characteristic functions, then the characteristic function of z will be given by
(9)QZ(ω)=Q1(ω)Q2(ω)⋯QN(ω).
So, the form of z can be simply obtained by an inverse Fourier transformation on the product of the characteristic functions.Define the characteristic functions forρk and ρk+1 as
(10)Qρk(ω)=∫-∞∞PDF(ρk)eiρkωdρk,Q-ρk+1(ω)=Qρk+1*(ω)=∫-∞∞PDF(ρk+1)e-iρ,k+1ωdρk+1.
Then making use of the above-stated method, the characteristic function Qδ(ω) of the required PDF(δk) is.
(11)Qδ,k(ω)=Qρk+1(ω)Qρ,k*(ω).
Then, required PDF(δk) can be obtained by simply using the convolution theorem:
(12)PDF(δk)=FT-1[Qρk+1(ω)Qρk*(ω)]=∫-∞∞PDF(ρk)PDF*(δk-ρk)dρk=∫-∞∞PDF(ρk)PDF(δk-ρk)dρk.
Hence, required PDF(ρk)  is the autocorrelation function of the positional disorder PDF(ρk). The integral is taken over the displacement of the metaatom from its mean position and limited by the finite values of the slot length. The strength of the method is the fact that no explicit assumption has been made regarding the form of  PDF(ρk) describing the positional disorder.The mathematical procedure has to ensure that the perturbation does not become so large that the meta-atoms overlap each other. In the analysis, the particles are assumed to be placed in average in the center of the slots of a length equal to the mean spacing period. The particles can randomly move within their own slot, and the extent of the displacement from the center of the slot is given byPDF(ρk). A consequence of such a restraint is that PDF(ρk) has to be restricted and normalized within this slot. Figure 5 shows how the autocorrelation function approaches a triangular function from its initial Gaussian form, as the position of the particle within the slot becomes completely random (i.e., PDF(ρk) takes a rectangular form). A simple algebraic form of the probability distribution function cannot be obtained due to this truncation. Hence, the following approach was adopted: the normalized versions of  PDF(δk) for the inter separation were obtained using numerical code, and they were subsequently used for numerical integration (according to (1)) for obtaining the effective material parameters.


## 3. Case of Randomly Positioned Dipoles

In this section, using the above mentioned principles, the effect of disorder in a chain of periodically placed dipoles is investigated. The geometry is given in Figure3. The bold arrow shows the direction of propagation of the electromagnetic wave.Figure 3
Geometry of propagation for randomly arranged dipole ensemble.The system can be mathematically modeled as follows. Considering the coupling dynamics between two equal adjacent oscillators, one can write the equation describing their dynamics as(13)∂2x1∂t2+γ∂x1∂t+ω02x1+σx2=qmEx,∂2x2∂t2+γ∂x2∂t+ω02x2+σx1=qmEx.
The term on the right side is the same for both oscillators, and as the same field impinges on both of them. By substituting the time ansatz xi(t)=xi(ω)exp(-iωt), the system can be easily solved for x1(ω) and x2(ω) as follows:
(14)[RσσR][x1(ω)x2(ω)]=qm[Ex(ω)Ex(ω)],R=ω02-ω-iγω.
Thus,
(15)x1(ω)=x2(ω)=q/mω02-ω-iγω+σEx(ω),σ=σ(z)=σ0a03z3.
Here, σ0 and a0 are the coupling constant and the distance between the oscillator. It is assumed that the interaction between the oscillators is the near field one between the dipoles that stipulates the inverse cubic distance dependence in the second equation in (15).The response of the system can thus be obtained by monitoring the susceptibility of the medium. The polarization of the system can be written as(16)Px(z,ω)=2q2ηm1ω02-ω-iγω+σ(z)Ex(ω),
so that the effective susceptibility is
(17)χx(z,ω)=2q2ηm1ω02-ω-iγω+σ(z).
To incorporate the effect of disorder, following (1), the averaged form of the above susceptibility can be obtained as
(18)〈χ(z0,ω)〉z=∫-∞∞PDF(δ,D)χ(z0,δ,ω)dδ
or
(19)〈χ(z0,ω)〉z=2q2ηm∫-∞∞PDF(δ,D)(z0+δ)3(ω02-ω-iγω)(z0+δ)3+σ0a03dδ,
where PDF(δ,D) is the inter separation PDF and D here quantizes the amount of disorder presented in the system.


## 4. Case of Randomly Positioned Quadrupoles

The extension of the above model to metamaterials (i.e., taking into account the magnetic response) firstly requires that the interaction between the adjacent meta-atoms is taken into consideration. The system is taken to be similar as the one shown in Figure3, but the dipoles are now replaced by quadrupoles (Figure 4). The long axis of the cut wires is oriented along the x-axis. The cut wires forming the quadrupole are separated along the y direction. The meta-atoms are arranged randomly (in terms of above described random positioning in the respective slots) along the z direction.Figure 4
Geometry of propagation for randomly arranged quadrupole ensemble.Relationship between the positional disorder function (a) and the interseparation probability distribution function (b). As the positional PDF (b) deviates from the Gaussian form for higher values of disorder (due to restrained excursion), the interseparation PDF approaches a triangular form.
(a)
(b)Assuming that a plane electromagnetic wave now propagates through the ensemble along they direction, while its electric vector is polarized along the x direction, the coupled dynamics of two meta-atoms can be modeled via four coupled oscillator equations:
(20)∂2x1∂t2+γ∂x1∂t+ω02x1+ax2+bx3+cx4=qmEx,1,∂2x2∂t2+γ∂x2∂t+ω02x2+ax1+bx3+cx4=qmEx,2,∂2x3∂t2+γ∂x3∂t+ω02x3+ax4+bx1+cx2=qmEx,3,∂2x4∂t2+γ∂x4∂t+ω02x4+ax3+bx2+cx1=qmEx,4,
where
(21)a=σ0,b=σ0y03z3,c=σ0y03(y02+z2)3/2,
and  σ0 is the value of the coupling constant measured for the interseparation y0. The magnitude of the coupling constant varies inversely as the cube of the distance, and so its value can be obtained for other interseparations—here z and the diagonal distance (y02+z2)1/2. The exponential phase factors in the right side take into account the retardation effect. It is clear that a change in the excitation conditions will affect the form of the right hand side of the above equations, while a change in the configuration of the meta-atoms can be accounted by a change in the form of the coupling coefficients. The procedure of determining the response of the medium then remains the same—one seeks to determine effective susceptibilities (corresponding to the symmetric and antisymmetric modes of oscillation), average them over all possible coupling configurations, and then use these values for ascertaining the effective material parameters.The first step is to find the solution of the above set of equations. They can be transferred to the Fourier domain by using ansatzxi(t)=xi(ω)exp(-iωt), i=1,2,3,4, so that the system can be rewritten in a matrix form:
(22)[RabcaRcbbcRacbaR][x1(ω)x2(ω)x3(ω)x4(ω)]=qm[Ex(ω)exp(iky1)Ex(ω)exp(-iky1)Ex(ω)exp(iky1)Ex(ω)exp(-iky1)],R=ω02-ω-iγω.
The modes of oscillation of interest are given by x1(t)±x2(ω). The system can be solved to obtain the values of x1(t)and x2(ω), and the modes of the system can be written as.
(23)x1(t)+x2(ω)=qExm2cos(ky1)ω02-ω-iγω+(a+b+c),x1(t)-x2(ω)=qExm2cos(ky1)ω02-ω-iγω-(a-b+c),
and hence one can define the effective susceptibilities:
(24)χ±(z0,δ,ω)=qExm2cos(ky1)ω02-ω-iγω±(a±b+c),
where the z0 and δ dependences are due to a,b, and c. Due to this form of definition, the functional forms of the polarization, quadrupole moment, and magnetizaton remain the same:
(25)P=2qηy1(2χ+(ω)cos(ky1)00)Ex(y,ω),Q=qηy1(02iχ-(ω)sin(ky1)02iχ-(ω)sin(ky1)00000)×Ex(y,ω),M=qηy1(002iχ-(ω)sin(ky1))Ex(y,ω).The effect of disorder can then be taken into account by carrying out the extra averaging integration (1):
(26)〈χ±(z0,ω,D)〉=∫-∞∞PDF(δ,D)χ±(z0,δ,ω)dδ,
or, more explicitly:(27)〈χ±(z0,ω,D)〉=∫-∞∞PDF(δ,D)1(ω02-ω-iγω)±σ0(1±(y03/(z0+δ)3)+(y03/(y02+(z0+δ)2)3/2))dδ,wherez0 is the mean period. The limits of the integration indicate that the autocorrelation procedure for PDF(δ,D) has been already carried out. This integral can be solved numerically for a given value of frequency.With the effective susceptibility as defined above, one may now consider a planar metamaterial, which is formed by the identical rows of the randomly positioned meta-atoms. The effect of randomness is taken into account by the averaging procedure, and hence the dispersion relation and the effective material parameters can be written in analogy to [16]. Adhering to the same conditions of geometry and excitation, the following expressions can be utilized:(28)ky2(ω)=ω2c21+A〈χ+(z0,ω,D)〉1+(ω2/c2)Ay12((1/2)〈χ+(z0,ω,D)〉-〈χ-(z0,ω,D)〉),εeff=1+A〈χ+(z0,ω,D)〉-Aky2y122(12〈χ+(z0,ω,D)〉-〈χ-(z0,ω,D)〉),μeff=11-(ω2/c2)Ay12〈χ-(z0,ω,D)〉.The previous expressions can be easily carried over to a numerical code to obtain the material parameters of interest. The following section presents the results and compares them with the experimental observations.


## 5. Method of Numerical Implementation

For convenience, the integrations and other expressions have been converted to their normalized versions. The frequencies are normalized with respect to the resonant frequencyω0 of the independent cut wire, while the distances are normalized with respect to the cut-wire spacing a0. Specifically, the susceptibilities for the case of dipoles and quadrupoles are(29)〈χ±(z0,ω)〉z=∫-∞∞PDF(δn,Dn)(1+δn)3(1-ωn2-i(ωn/Q))(1+δn)+σ0/ω02dδn,〈χ±(z0,ω)〉z=∫-∞∞PDF(δn,Dn)(1-ωn2-i(ωn/Q))±σ0(1±1/((z0)n+δn)3+1/(1+((z0)n+δn)2)3/2)dδn,where(30)δn=δa0,Dn=Da0,ωn=ωω0PDF(ρk,D) is assumed to be Gaussian:
(31)PDF(ρk,D)=12Dexp(-ρk22D2),-a02<ρk<a02,PDF(ρk,D)=0,-a02<ρk,ρk>a02.
The excursion of the dipoles around their mean positions had to be limited within the interval [-a0/2,a0/2]; recalculation PDF(ρk,D) into PDF(δk,D) is given by (12). The integrals cannot be taken analytically and was done using the mathematical software MATLAB. Truncation of the positional PDF was achieved by coding. To obtain the autocorrelation of the PDF, a standard subroutine was used. The results of the operations are shown in Figure 5.All the constants used in the analysis were taken from [16]. The spacing between the cut wires a0 was taken to be 65 nm, and the resonant frequency of an isolated cut wire was taken as ω0=1.39×1015 rad s−1. The damping coefficient was taken to be γ=9.42×1013 rad s−1. The mean periodic spacing z0=1.8 (the mean spacing between the meta-atoms was taken to be 1.8 times the cut wire spacing y0).To verify the correct functioning of the code, the results for a very small disorder were compared with the result for a perfectly ordered system (with neighboring meta-atoms interacting with each other); see Figure6.Verification of the code—typical values from [16] were used in the computer code written for the calculation of the effective material parameters for metamaterials with positional disorder having a very small amount of disorder (D=0.01). The results obtained match with those in [16]—this is expected as the nature of coupling considered in the present theory should have negligible influence upon the material parameters for very large spatial periods (zn>3).
(a)
(b)


## 6. Results

The results of the analysis for dipoles are presented in Figure7, and the results of the analysis for quadrupoles are presented in Figures 8 and 9.Effective material parameter curves for dipole ensembles exhibiting positional disorder. The effective permittivity and permeability curves for disordered dipole ensembles are presented for different values of disorder. The first column pertains to values obtained for a mean period ofzn=1.2, while the second column pertains to those obtained for a mean period of zn=1.8. For the respective periodicities: (a) and (b) the positional disorder function; (c) and (d) the respective interseparation PDFs; (e) and (f) scaled real part of the permittivity; (g) and (h) scaled imaginary parts of the permittivity. Clearly, increase in disorder brings about a fall in the maximums of the response of the system.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Dispersion and effective material parameter curves for quadrupole ensemble withzn=1.2: (a) positional disorder function PDF; (b) interseparation PDF; (c) and (d) real and imaginary part of k-vector (in normalized units of cut-wire separation distance a0/2); (e) and (f) real and imaginary parts of effective permittivity; (g) and (h) real and imaginary parts of effective permeability.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Dispersion and effective material parameter curves for quadrupole ensemble withzn=1.8: (a) positional disorder function; (b) interseparation PDF; (c) and (d) real and imaginary part of k-vector (normalized with the cut-wire separation distance a0/2) (e) and (f) real and imaginary parts of effective permittivity; (g) and (h) real and imaginary parts of effective permeability.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)The analysis was carried out for two values of the spacing periodzn, namely, zn=1.2 (Figure 8) and zn=1.8 (Figure 9). The positional PDF was taken for the four different values of the standard deviation D, and consequently the inter separation PDF(δn,D) was obtained using numerical coding in MATLAB. The effective susceptibilities were calculated by numerical implementation of the integration (29) and (31), and then the effective material parameters (28) were calculated.The following features are clearly noted.(i)
For the disordered dipole ensemble, the fall in the permittivity with increasing disorder is clearly visible (Figures7(e), 7(f), 7(g), and 7(h)). The Im(ε) curve is symmetric for very small values of the variance (here, D=0.01). But as the disorder increases, the peak shifts towards lower frequencies, and the curves becomes broadened and asymmetric. The reason for these observed effects can be explained as follows. A disordered system can be thought to be made up of several different periodic systems. The resonant frequency for each periodic ensemble depends inversely on its spatial period. If the response of the disordered system is approximated by the sum of the responses of its constituent periodic systems, it becomes evident that the final curve will develop a tail approaching the blue end of the spectrum. The asymmetry can thus be attributed to an inverse power relationship between resonance frequency and interseparation. The effect of broadening is a consequence of particle conservation. On the other hand, the lowest frequency/largest wavelength of response is not a function of the periodicity, but is actually limited by the eigenfrequency of the independent oscillator. In fact, the resonance frequency approaches the eigenfrequency for a periodic assembly of dipoles, when the spatial period becomes large. Hence, as the disorder in the system increases, the curves become broadened and asymmetric, and the peak response shifts towards the eigenfrequency of the independent oscillator.(ii)
In the case of the quadrupole ensemble, a decrease in the value of the electric permittivity is observed asD is increased. This is in agreement with the experimental results. However, there is also a decrease in the value of the magnetic permeability. This decrease is more pronounced for zn=1.2 as compared to zn=1.8. This is an unexpected result, as the magnetic response should remain almost unaffected. The reason for this discrepancy could lie in the simple form of the probabilistic model chosen to describe the randomness.(iii)
Generally speaking, the final expressions for the permittivity and the permeability were derived under several approximations, associated with (1). The observed discrepancy could also be attributed to these approximations. Above all, the fundamental limitations of the multipole theory itself could affect the final results as well. These possibilities have to be investigated further.In the light of the above arguments, it is concluded that as the observed positions of the resonances and the relative magnitudes of the parameters are within the limits of approximation, the analysis is valid and can be used to roughly predict the properties of metamaterials with incorporated randomness.


## 7. Other Forms of Disorder

In the preceding analysis, the effect of positional disorder (arising due to aperiodicity) on the averaged material parameters was considered. In a random metamaterial other forms of disorder can exist as well. A particular case of interest is positional disorder along the cut-wire axis; see Figure10.Figure 10
Disorder along the cut-wire axis direction—the figure shows a one-dimensional disorder arrangement of meta-atoms. The extent of disorder can be quantified in terms of the angle, the total range of variation being limited to(-π/2,π/2).If this form of positional disorder is taken into consideration along with the aperiodicity, the model would then be a step closer to emulate a true self-organized random metamaterial [46]. In the multipole model, the individual cut wires are replaced by dipoles. In case the quadrupoles are disarrayed, the coupling between them will also be a function of their relative angular positioning. This angular dependence can be introduced into the coupling terms of the dynamic equations. More specifically, the coupling constants b and c in the differential equations will include the angular dependence. All other mathematics remaining the same, the averaging procedure can now be carried out between angles (-π/2,π/2).The curves in Figures12 and 13 summarize the results obtained by using the multipole approach.Two forms of the distribution function were used in the analysis (Figure11). The extent of disorder is correlated to the relative angular position of the dipoles. The first form of the angular PDF distribution function used was a rectangular function (see Figure 11(a)). The second form used was a Gaussian distribution (Figure 11(b)), the random variable being the relative angular position. The function is centered about 0 degrees, the extent of disorder being quantified by the standard deviation D; b and c are multiplied by the term cos(θ) to incorporate the angular dependence. Clearly then, when θ=π/2, there is no interaction between the cut wires. The Riemannian integration is limited between the values (-π/2,π/2). The constants were again taken from the original reference [22], and the mean periodicity was taken to be zn=1.8. The results (Figures 12 and 13) show that both the effective permittivity and permeability are clearly affected by the angular disorder. In a similar fashion, in-plane and out-of-plane skew disorders of meta-atoms can also be accounted by the model, by making appropriate changes to the coupling terms.(a) Rectangular and (b) Gaussian forms of distribution function used for governing the positional disorder of the metamaterials along the lateral direction. The positional disorder is expressed in terms of the relative angle between two neighboring meta-atoms. The effective material parameters of the ensemble are derived and presented in Figure12 for three different values of disorder.
(a)
(b)Effective material parameters for metamaterials exhibiting positional disorder along the lateral direction (y direction) governed by the rectangular distribution function (Figure 11(a)). The first column gives the material parameters for a metamaterial ensemble having a mean period of zn=1.2, while the second column is for metamaterials having a mean period of  zn=1.8.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Effective material parameters for metamaterials exhibiting positional disorder along the lateral direction (x-direction) governed by the Gaussian distribution function (Figure 11(b)). The first column gives the material parameters for a metamaterial ensemble having a mean period of  zn=1.2, while the second column is for metamaterials having a mean period of  zn=1.8.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)


## 8. Conclusions

In extending the multipole approach [16] to the case of random metamaterials, the effect of spatial distribution of the meta-atoms was taken into account by considering the near field coupling between neighboring meta-atoms. The effective susceptibility was expressed as a function of the inter-separation between meta-atoms, the ensemble averaged susceptibility was then obtained as an expectation value, weighed by the probability distribution function of all possible inter-separations. In the present work, the disorder was considered only along one direction (in-plane, perpendicular to cut-wire long axis), assuming that adjacent rows of meta-atoms do not interact with each other. Results obtained by the numerical implementation of the equations indeed confirm the experimental finding that increasing disorder has a more pronounced effect on the effective electrical permittivity than on the effective magnetic permeability. For smaller periodicities, however, the electrical permittivity and magnetic permeability are affected equally. This conflicting result could be caused by a coupling of the quadrupole moments of neighboring meta-atoms. This has not been explicitly considered in the present version of the model. Also, other factors such as the effect of incident polarization, or the coupling between adjacent rows, have not been considered in the present theory. The understanding gained from the study of this simple case can now be used to account for the above specific and more involved cases.



---

*Source: 161402-2012-11-25.xml*