Table Of Content

NONLINEAR MAGNETOINDUCTIVE TRANSMISSION LINES NIKOS LAZARIDES†, VASSILIS PALTOGLOU‡ and G. P. TSIRONIS⋆ Department of Physics, University of Crete, and FORTH P.O. Box 2208, 710 03, Heraklion, Crete, Greece †[email protected] ‡[email protected] ⋆ [email protected] Power transmission in one-dimensional nonlinear magnetic metamaterials driven at one end is investigatednumerically andanalytically inawidefrequencyrange. Thenonlinearmagneticmeta- 1 materials are composed of varactor-loaded split-ring resonators which are coupled magnetically 1 through their mutual inductances, forming thus a magnetoiductive transmisison line. In the lin- 0 ear limit, significant power transmission along the array only appears for frequencies inside the 2 linear magnetoinductive wave band. We present analytical, closed form solutions for the magne- n toinductive waves transmitting the power in this regime, and their discrete frequency dispesion. a When nonlinearity is important, more frequency bands with singificant power transmission along J the array may appear. In the equivalent circuit picture, the nonlinear magnetoiductive transmis- 3 sion line driven at one end by a relatively weak electromotive force, can be modeled by coupled 1 resistive-inductive-capacitive(RLC)circuitswithvoltage-dependentcapacitance. Extendednumer- ical simulations reveal that power transmisison along the array is also possible in other than the ] linearfrequencybands,whicharelocatedclosetothenonlinearresonancesofasinglenonlinearRLC h circuit. Moreover, the effectiveness of power transmission for driving frequencies in the nonlinear p bands is comparable to that in the linear band. Power transmission in the nonlinear bands occurs - s throughthelinear modesofthesystem,andit isclosely related totheinstability of amodethatis s localized at thedriven site. a l c PACSnumbers: . Keywords: nonlinear split-ring resonators, magnetic metamaterials, driven linear waves, nonlinear power s transmission c i s y I. INTRODUCTION materials. Magnetism at such high frequencies is partic- h ularly important for the implementation of devices such p as tunable lenses, adaptive mirrors, etc., since there are [ The metamaterial concept usually refers to a periodic only a few natural materials that respond magnetically arrangementof artificially structured elements, designed 1 above microwaves. Moreover, magnetism in those natu- v to achieveadvantageousand/orunusualelectromagnetic ralmaterialsisusuallyweakandwithinnarrowfrequency 9 propertiescomparedtothoseofnaturallyoccuringmate- bands, limiting thus their use in possible Terahertz de- 1 rials [42]. For example, the first metamaterial ever real- vices. The most common realization of such a magnetic 5 ized [40], exhibited negativerefractionindex in a narrow 2 metamaterial iscomposedofperiodicallyarrangedSRRs frequency band around 5 Gigahertz. The negative re- . in a one- ortwo-dimentionallattices. In its simplest ver- 1 fraction index property requires that both the dielectric sion an SRR is just a highly conducting [5, 13, 24] or 0 permittivity and the magnetic permeability are simul- 1 superconducting [8] ring with a slit. The SRR structure taneously negative. The most widely used elements for 1 wasoriginallyproposedbyPendry[30],whotheoretically constructinganegativerefractiveindexmetamaterialare : predicted that a periodic SRR-based structure could ex- v the electrically small resonant’particles’called split-ring hibit a negative magnetic permeability frequency band i resonators (SRRs)andmetallicwires. Whilethearrayof X when irradiated by an alternating electromagnetic field wiresgivesa negativepermittivity below the plasma fre- r ofappropriatepolarization. Pendryalsorealizedthatthe a quency, the arrayof the SRRs gives a negative magnetic SRR structure has considerable potential for enchancing permeability above their resonance frequency. Thus, the nonlinear phenomena. When the SRR dimensions are refraction index can be negative in a narrow frequency much smaller than the wavelength of the incident field, band. Sincetheirdiscovery,metamaterialshavebeenthe the SRR canbe regardedasanelectric circuitconsisting target of intensive study. Although the developement of ofinductive,capacitiveandohmicelementsconnectedin metamaterials at microwave frequencies has progressed series. The ring forms the inductance L, while the slit tothepointwherescientistsandengineersarenowpursu- can be considered as a capacitor of capacitance C. The ing applications, researchon metamaterials that operate ohmic element R may represent all possible losses of the athigherfrequenciesis stillinanearlystage[25,26,43]. structure. Thereisawidesubclassofmetamaterials,thatexhibit significant magnetic properties and even negative mag- Thenextstepwasnaturallytobuiltnonlinearityinthe netic permeability at Terahertz and optical frequencies SRR structure. Nonlinearity implies real-time tunability [13, 25, 51], even though they are made of non-magnetic and multistability that is certainly a desired property of 2 possiblefuture devices. TheSRRs canbecomenonlinear toinductive wave modes, slightly modified by the pres- either by the insertion of a strongly nonlinear dielectric enceofnonlinearity. However,forfrequenciesinthepass- [9] or a nonlinear electronic component (i.e., a varactor) bandsresultingfromnonlinearity,the poweris transmit- [33,34,50]intotheirslits. Bothwayswouldleadtoeffec- ted along the array from nonlinear, breather-like exci- tivelyfielddependentmagneticpermeabilityforanarray tations. In the next Section we derive the model equa- of nonlinear SRRs. Recently, the dynamic tunability by tions for a varactor-loaded SRR-based magnetic meta- anexternalfieldfora two-dimensionalarrayofvaractor- material, where the nonlinearity has the form a polyno- loadedSRRswasdemonstratedexperimentally[37]. The mial expansion on the main variable, like that investi- SRRs in such an array are coupled through magnetic gated in [50]. We also present bifurcation diagrams for and/orelectricdipole-dipoleinteractions,whosestrength the shortestpossibleend-drivenarray,i.e.,anarraywith dependsontherelativeorientationoftheSRRsandtheir only two varactor-loaded SRRs, from which only one is slits in the array [6, 10, 31, 45]. directly driven. In Section III we obtain analytical solutions for the linear magnetoinductive waves of the end- In the equivalent circuit picture, the dynamics of a driven,finite SRRarraywhenlossesareneglected,along periodic nonlinear SRR array irradiated by an electro- with their discrete dispersion relation. In Section IV we magnetic field, can be described by a set of coupled and presentnumericalresultsforpowertransmissionanddis- driven ordinary differential equations with on-site non- cuss their dependence on the model parameters for a linearity. The discreteness, which is inherent in those short, end-driven varactor-loaded SRR array with and SRR-basedmagneticmetamaterials,alongwiththe non- without an absorbing boundary at the non-driven end. linearity and the weak coupling between their elements, We finish with concluding remarks in Section V. allow for the generation of nonlinear excitations in the form of intrinsic localized modes or discrete breathers [41]. These nonlinear modes appear generically in dis- II. VARACTOR-LOADED SRR ARRAY MODEL crete and nonlinear extended systems. Recent theoreti- calworkinone-andtwo-dimensionalnonlinearSRRlat- ticeshavedemonstratedtheexistenceandthestabilityof Consider a ring-shaped or a squared split-ring res- discretebreathersbothinenergy-conservinganddissipa- onator with a hyperabrupt tuning varactor mounted tive systems [3, 4, 16, 18]. It has also demonstratedthat onto its slit. The varactor type could be selected breathers may be formed spontaneously through modu- from a large variety of available varactors, according lational instability in binary nonlinear magnetic meta- to the needs of the experiment. In recent experimen- materials[20,22,29]. Moreover,domain-wallexcitations tal works on varactor-loaded SRR-based metamateri- [35] and envelope solitons [1, 2, 14, 48] may as well be als [11, 15, 32, 50], the varactor selected was a Sky- excited in those systems, which seem to be stable even works SMV1231-079, whose voltage-dependent capaci- in the presence of noise [49]. tanceC(UD)(UD isthevoltageacrossthediode)isgiven by The SRR-based magnetic metamaterials support a new kind of electromagnetic waves, the magnetoinduc- U −M D tive waves, which exhibit phonon-like dispersion curves, C(UD)=C0 1 , (1) − U andthey cantransferenergyalongthe array[39,46,47]. (cid:18) p(cid:19) It is thus possible to fabricate a contact-free data and where C0 is the DC rest capacitance, Up is the intrinsic power transfer device, a magnetoinductive transmission potential,andM isaparameter. Thevaluesofthosepa- line, which make use of the unique properties of the rametersareprovidedbythemanufacturerSPICEmodel magnetic metamaterial structure, and may function as for that varactor to be a frequency-selective communication channel for devices M =0.8, C =2.2 pF, U =1.5 V. (2) via their magneto-inductive wave modes [44]. In the 0 p presentworkweinvestigatethepowertransmissionalong From Eq. (1) we can determine the voltage-dependence a one-dimensional varactor-loaded SRR-based magnetic of the normalized charge metamaterial, in the two possible geometries, which is Q driven at one end by a sinusoidal power sourse. In the q = D , (3) C U linearlimit,energytranferalongthearrayoccursonlyfor 0 p driving frequencies in the linear magnetoinductive wave with Q being the charge in the diode, as a function of D band. However,thenonlinearitycouldgeneratemorefre- the normalized diode voltage u=U /U as D p quency bands where efficient power transmission along 1 the array is possible. Power transmission in chains of q = 1 (1 u)−M+1 . (4) coupled anharmonic oscillatorsfor driving frequencies in M +1 − − − h i the band gap of the linear spectrum has been recently Assuming that U < U , the earlier equation can be D p investigated,andthateffectis referredto asself-induced solved for the diode voltage as a function of the charge, transparency [28] or supratransmission [7]. For frequen- i.e., u=u(q) cies inside the linear magnetoinductive wave band the 1 power is transmitted along the array from the magne- u=1 [1 q( M +1)](−M+1) . (5) − − − 3 When x q( M+1) <<1,thevoltageuintheearlier 5 | |≡| − | equation can be expanded in a Taylor series inpowers of q around zero. Thus, neglecting term of order (q4) or 4 O higher, we get 3 u=q+αq2+βq3, (6) 0 i where, by using Eq. (2), 2 M 1 1 α= = 0.4, β = M(2M 1)=+0.08. (7) − 2 − 6 − The expression Eq. (6) can then be used in the volt- 0 0 1 2 age equation, obtained from Kirchhoff’s voltage law, Ω whichrepresentsthevaractor-loadedSRRasaneffective resistive-inductive-capacitive(RLC)circuitwithvoltage- FIG.1: Thecurrentamplitudei0asafunctionofthenormal- dependent capacitance ized frequency Ω for a varactor-loaded SRR modeled by Eq. (13). Red curve: frequency increases; black curve: frequency d2Q dQ D D decreases. L +R +U u= (t), (8) dt2 dt p E where driven by a sufficiently strong external field. A typical Q Q 2 Q 3 current amplitude i0 - frequency Ω curve for the equiva- U u U =U D +α D +β D , lentcircuitmodelEq. (13),likethatshowninFig. 1,ex- p D p ≡ "C0Up (cid:18)C0Up(cid:19) (cid:18)C0Up(cid:19) # hibitsallthosecharacteristics. Noticethehysteresisloop (9) at normalized frequency Ω around 1.1, close to the reso- , (t) is the induced electromotive force of amplitude 0 nancefrequencyofthelinearsystem,andtheverystrong E E and frequency ω resulting from an applied electromag- resonance at its second harmonic, i.e., at Ω 2, where netic field with the same frequency, of the form a smaller hysteresis loop also appears. Furth≃ermore, we observe a subharmonic resonance at Ω 0.48, as well as E(t)=E0 sin(ωt), (10) weaker resonances at Ω 0.63, 0.33 an≃d Ω 0.25. For ∼ ∼ frequency intervals lying between the boundaries of the where L is the inductance of the ring, R is the sum of hysteresis loops there are two different states of the os- the Ohmic resistance of the SRR ring and the series re- cillator, with low and high current amplitude, which are sistance of the diode, and t the temporal variable. The simultaneously stable. eigenfrequency ω of the circuit described by Eq. (8) is 0 given, in the linear limit without losses and driving, by H D 1 E ω = . (11) w 0 √LC0 h l d That approximate model of the varactor-loaded SRR is H foundtobe adequateforvaractorvoltagesnotexceeding E 0.5V (u 1/3innormalizedunits),whenthedissipative ≤ current of the varactor can be neglected [50]. Using Eq. (3) and the relations FIG. 2: Schematic of an SRR array in the planar geometry ω (upper)and in the axial geometry (lower). 0 τ =tω , Ω= , ε = E , (12) 0 0 ω U 0 p Consider a one-dimensional periodic array of identical Eq. (8) can be written in the normalized form varactor-loadedSRRs,whicharecoupledmagneticallyto their first neighbors through their mutual inductances. d2q dq +q+αq2+βq3+γ The array can be formed in two possible geometries, ac- dτ2 dτ cordingtotherelativeorientationofitselementaryunits =ε0sin(Ωτ), (13) (i.e., the SRRs) in the array. As can be seen in Fig. 2, the SRRs can be arranged either in the planar geome- where try, where all of them are lying on the same plane, or in the axial geometry, where the axes of all the SRRs are γ =R C /L (14) 0 lying on the same line. The orientation of the incident is the normalized loss coeffipcient. The varactor-loaded electromagneticfield,whichisalsoshowninthefigure,is SRR is a nonlinear oscillator, which exhibits multista- suchthatthe magnetic componentcanexcite anelectro- bility, hystereticeffects,andsecondaryresonances,when motiveforceineachSRR. Belowweconsiderarraysthat 4 aredrivenonlyatoneend(letus saytheleftend). Thus 20 (a) only one SRR receives energy directly from the driver. 15 Hduoewetovetrh,ethienteenrearcgtyiocnanbebtweeterannistmsietlteemdeanlotsn.gItnhceluadrrinagy Etot10 thecouplingofeachSRRtoitsnearestneighbors,weget 5 for the charge Q of the n th varactor the dynamic D,n − 0 equations 1 (b) dI dI dI n n−1 n+1 L +RIn+Upun+M +M = δn,1(,15) 2 dt dt dt E f 0.5 , 1 where n = 1,...,N, with N being the total number of f SRRs in the array, I (t)=dQ /dt is the currentflow- n D,n 0 ing in the n th SRR, U u is the approximate voltage p n across the n−th varactor, and M (not to be confused 1 (c) − with the varactor parameter) is the mutual inductance 2 betweentwoneighboringSRRs. Usingtherelationsfrom f 0.5 , Eq. (12), thecouplednonlinearEqs. (15)canbewritten 1 f in normalized form as 0 d2 1 2 dτ2 (λqn−1+qn+λqn+1)+qn+αqn2 +βqn3 = Ω dq n γ +ε sin(Ωτ)δ , (16) − dτ 0 n,1 FIG.3: (a) Energy Etot bifurcation diagram forvaryingfre- quency Ω for an asymetrically driven varactor-loaded SRR whereλ=M/Listhecouplingcoefficientbetweenneigh- dimer with γ = 0.001, ε0 = 0.3, and λ = +0.01. (b) Energy boring SRRs. The function δn,1 is unity for n=1, while fractions f1 and f2 of thefirst (driven) and the second SRR, itiszeroforanyothern. Thecouplingparameterλmay respectively,ofthedimerinthelow energystate. (c)Energy assume bothpositive and negativevalues,corresponding fractions f1 and f2 of thefirst (driven) and the second SRR, toplanarandaxialgeometryofthevaractor-loadedSRR respectively, of thedimer in thehigh energy state. array,respectively. For a particular array, the values of the coupling and loss coefficients λ and γ, respectively, can be estimated (19)thatRSRR 0.13Ω. Thediode seriesresistanceRs ≃ from its geometrical and material parameters. Specifi- of the varactor used in Ref. [50] is around 2 Ω; that ≃ cally, the mutual inductance M, the ring inductance L, is, most of the Ohmic resistance of the varactor-loaded and the Ohmic resistance RSRR of the ring can be cal- SRR,R=RSRR+Rs =0.13+2.00=2.13Ωcomesfrom culated from the following expressions the varactor. By substituting R 2.13 Ω in the second ≃ of Eq. (14), along with L=14 nH and C =2.2 pF, we 0 πa2a2 get that γ = 0.027. The resonance frequency of a single M = µ , (17) 0 4D3 varactor-loaded SRR is f0 (2π√LC0)−1 0.9 GHz. ≃ ≃ 16a In the numerical calculations below we used the esti- L = µ a ln 1.75 , (18) 0 mated value of λ that indicates weakly coupled SRRs. h − (cid:20) (cid:18) (cid:19) (cid:21) However, we use a smaller value of γ, which could be 2a R = , (19) achieved in practice by using a metalic wire with higher SRR σhδ conductivity and an ultra-low resistance varactor with where µ is the magnetic permeability in vacuum, a is similar capacitance-voltagerelation. 0 the average radius of each SRR, h the diameter of the The shortest end-driven varactor-loadedSRR array is wire of each SRR, D is the center to center distance be- thatcomprisedoftwoelements. Sinceonlyoneofthemis tweenneighboringSRRs,andσ andδ istheconductivity drivenbytheexternalfield,werefertothissystemasthe and skin depth of the SRR wire. Eqs. (17)–(19) refer to asymmetrically driven nonlinear dimer. The increase of circular SRRs with circular cross-section. However,they degrees of freedom leads to interesting and complicated can also be used for circular SRRs with square cross- behavior, that is marked by the appearance of chaos. In section by defining an equivalent diameter of the SRR Fig. 3,thebifurcationdiagramofthetotalenergyE of tot wire h = 4dw/π, with d and w being the depth and the dimer with varyingnormalizedfrequency Ω is shown width, respectively, of the SRR wire. For a 3.7 mm, for positive coupling λ. For negative coupling λ with p ≃ D 11 mm, h 0.5 mm, we get from Eqs. (17) and equal strength we get very similar results that share a ≃ ≃ (18), respectively, that M 0.14 nH and L 14 nH, number of characteristic features. First of all, multista- and consequently λ 0.01≃. For δ 2.2 ≃10−3 mm bility is again present, for normalized frequencies in an and σ 0.5 106 (Ω≃cm)−1 (for an SR≃R mad×e of copper interval around Ω 0.9. There, two simultaneously sta- ≃ × ∼ wire at 1 GHz at room temperature) we get from Eq. bleenergystates,withhighandlowenergy,co-exist,asit ∼ 5 isshowninFig. 3a. Theenergydifferencebetweenthose the structure, i.e., states increases condiderably with increasing frequency. Moreover,low energy state becomes chaotic for frequen- q0 =qN+1 =0. (22) cies around Ω 1.9. It is also interesting to see how the ∼ By substituting q =Q sin(Ωτ) into Eq. (20), the sta- total energy is divided between the two elements. The n n tionary equations can be written as energy fractions f = E /E and f = E /E , with 1 1 tot 2 2 tot E and E being the energies of the first (driven) and 1 2 sQ +Q +sQn+1=κ, (23) n−1 n second SRR of the dimer, are shown both for the high and the low energy states in Fig. 3b and 3c as a func- where tion of the normalized frequency Ω. When the dimer is in the low energy state, most of the energy is concen- λΩ2 ε0 s= , κ= δ , (24) trated in driven SRR whose energy fraction is close to −1 Ω2 1 Ω2 n,1 − − unity (correspondingly, the energy fraction of the other or, in matrix form SRR is close to zero). However, for frequencies around the resonances, energy can be transfered easily from the Q=κSˆ−1E , (25) 1 one SRR to the other, and then the energy fractions of the SRRs in the dimer attain comparable values. In the where Q and E are N dimensional vectors with com- 1 chaotic regions, the energy is tranfered irregularly from ponets Q and δ , resp−ectively, and Sˆ−1 is the inverse i i,1 one SRR to the other and the energy fractions span the of the N N coupling matrix Sˆ. The latter is a real, whole interval between zero and one. In the high energy symmetric×tridiagonalmatrixthathasdiagonalelements state (Fig. 3c), the energy fractions are almost equal for equal to unity, while all the other non-zero elements are frequencies above the linear resonance frequency, Ω 1. equal to s. The need to find the inverse of tridiagonal ∼ Belowthatfrequency,theenergyiseitherconcentratedin matrices like Sˆ arises in many scientific and engineering the driven SRR (for frequencies far from resonances), or applications. Recently, Huang and McColl [12] related the values of the energy fractions are comparable (close the inversion of a general trigiagonal matrix to second to resonances. order linear recurrences,and they provided a set of very simpleanalyticalformulaefortheelementsoftheinverse matrix. Thoseformulaeleadimmediatelytoclosedforms III. DRIVEN LINEAR SOLUTIONS for certain trigiagonal matrices, like Sˆ [21]. The components of the Q vector can be written as For a weakly driven array, the nonlinear terms which are proportional to α and β are not so important and Q =κ Sˆ−1 , (26) i they can be neglected. Then, the system of Eqs. (16) i1 (cid:16) (cid:17) reduces to a linear one where Sˆ−1 is the (i1) element of the inverse cou- d2 (λq q +λq )+q pling m(cid:16)atrix(cid:17)iSˆ1−1, whose e−xplicit form is given in [21]. dτ2 n−1− n n+1 n Then, the solution of the linear system Eq. (20) and dq n (22) for any driving frequency and finite N is given by = γ +ε sin(Ωτ)δ . (20) 0 n,1 − dτ sin[(N n+1)θ′] The dispersion relation of the linear system can be ob- qn(τ)=κµ sin[(N−+1)θ′] sin(Ωτ), tained by the substitution q = Acos(κn Ωτ) in the n − 1 absence of losses and applied field, i.e., for ε0 = 0 and θ′ =cos−1 , (27) γ =0, which gives the linear dispersion of magnetoiduc- 2s (cid:18) | |(cid:19) tive waves for s>+1/2 and s< 1/2, and − 1 Ω = , (21) κ sinh[(N n+1)θ] 1+2λ cos(κ) q (τ)=κµ − sin(Ωτ), n sinh[(N +1)θ] p where κ is the (normalized) wavevector ( π κ π). 1+ 1 (2s)2 Eq. (21) defines a frequency band of wi−dth≤∆Ω ≤ 2λ θ =ln − , (28) thatisboundedbyaminimumandamaximumfrequ≃ency p2|s| Ωmin =1/√1+2λandΩmax =1/√1 2λ,respectively. for 1/2<s<+1/2, where However,itis alsopossibleto calcula−te the dispersion, − expressedas a series of resonantfrequencies, and the ex- 1 s n−1 µ= | | . (29) act form of the linear magnetoinductive modes for the s − s finite system in the absence of losses (γ = 0). For finite | |(cid:18) (cid:19) systems, Eq. (20) should be implemented with free-end Note that the Q s are uniquely determined by the pa- i boundary conditions, to account for the termination of rameters of the system. 6 20 (a) n 10 150 Q 10 0 8 100 40 (b) 50 6 Qn 0 n 0 4 -40 -50 40 (c) 2 -100 Qn 0 00.99 1 1.01 -150 -40 Ω 40 (d) n Q 0 FIG.4: DensityplotsoftheQn’sinastationarystate,onthe sitenumbern-normalizedfrequencyΩplanecalculatedfrom -40 the analytical expressions Eqs. (27) and (28) for λ = +0.01 20 andε0 =0.3. Theverticalwhitebandsindicatethepositions n 0 of the resonances. Q (e) -20 From the analytical solution in the linear magentoin- 2 4 6 8 10 n ductive wave band, Eq. (27), which corresponds to either s > +1/2 or s < 1/2, we can infer the resonant − frequenciesfromzeroingthedenominator,i.e.,bysetting FIG. 5: Typical profiles calculated analytically from Eqs. sin[(N +1)θ′] = 0. Then we see that the solution has a (27) and (28) for λ = +0.01, ε0 = 0.3, and five different frequencies: (a) Ω=0.989; (b) Ω=0.995;(c) Ω=1.000; (d) resonance for Ω = 1.005; (e) Ω = 1.011. The frequencies in (a) and (e) lie 1 outside thelinear magnetoinductivewave band. s s = , (30) m ≡ 2cos mπ (N+1) h i see in Fig. 4 that the solutions close to the lower bound where m is an integer (m = 1,...,N). The resonant fre- (Ω 0.99)ofthe linearmagnetoinductivewaveband min qquencies are obtained by solving the first of Eq. (24) exhibit≃very high amplitudes. For negative coupling co- withrespecttoΩ,andsubstitutingtheresonantvaluesof efficientwithequalstrengththehighamplitudesolutions sdi≡spesrmsifornomreElaqti.o(n30o)f.lTinheeanr wmeaggentettohiendduiscctrievteewdyanveasmiinc alipnpeaerarmcalogsneettoointdhuectuipvepewrabvoeubnadnd(.Ωmax ≃ 1.01) of the the varactor-loaded SRR array, the discrete analogue of InFig. 5theprofilesoftheanalyticalsolutionsforfive Eq. (21), as different frequencies, both inside and outside the linear magnetoinductive wave band are shown for positive λ. 1 Ω = , (31) The frequencies in Figs. 5a and 5e lie outside the lin- m 1 2λ cos mπ ear band, so that the Qn’s decrease exponentially with − N+1 increasing site number n. r (cid:16) (cid:17) wheremisthemodenumber(m=1,...,N). InFig. 4,a densityplotoftheQ ’sonthen Ωplane,calculatedan- n IV. POWER TRANSMISSION ALONG THE − alytically from Eqs. (27) and (28), is shown for positive ARRAY. coupling coefficient λ. In the frequency interval shown, theQ ’shaveratherlargevaluesbecause,duetothenar- n In the following we focus on the power transmission row band, all frequencies are very close to a resonance. along a varactor-loaded SRR array with N = 10 ele- Moreover, the analytical solutions presented above have ments, which is driven at the left end (n = 1). The been obtained in the absence of losses, so that there is average power dissipated in the n th SRR is defined as nothing to prevent the Q ’s from increasing indefinitely − n closetoaresonance(andtogotoinfinityexactlyatreso- P =R<I2(τ)> , (32) n n τ0 nance). The white verticalregionsin these figures corre- spond to frequency intervals around a resonance, where where <> denotes time-average over τ time-units. τ0 0 the Q ’s have very large values. Thus, those regions in- The power density P can be normalized to P = I2R, n n 0 0 dicate the frequencies of the resonances. We can also where I = ω C U and R the value of resistance that 0 0 0 p 7 results from the value of the dissipation coefficient used in the simulations. For γ =0.001 we have that R 0.06 ≃ andP0 28µW. Itis convenientto expressPn indBm, 10 0 ≃ according to the relation 8 -50 P (in W) P (dBm)=20log <i2 > 0 n (cid:18) n τ0 1 mW (cid:19) n 6 -100 ≃4.34ln(0.028<i2n >τ0), (33) -150 4 where i is the normalized current in the n th SRR. n − -200 2 -250 10 50 0 0.5 1 1.5 2 -300 f (GHz) 0 8 -50 n 6 FIG. 7: Average power density plot (in dBm) on the n−f -100 plane, for λ = +0.01, γ = 0.001, α = −0.4, β = 0.08, and 4 -150 ǫ0 =0.3. -200 2 -250 thelosscoefficienthasbeenincreasesbyanorderofmag- nitude (from 0.001 to 0.01. In the density plots shown 0 0.5 1 1.5 2 -300 above,significant power transmissionis indicated by red f (GHz) color. In both Figs. 7 and 8, nonlinear power transmission occurs at f 0.3 GHz, f 0.45 GHz, and ∼ ∼ f 1.8 GHz. The latter nonlinear band, around the FIG. 6: Average power density plot (in dBm) on the n−f ∼ secondharmonicofthe resonancefrequencyofthe linear plane, for λ = +0.01, γ = 0.001, and ǫ0 = 0.3, in the linear system, seems to be the wider of all. regime (α=0, β=0). The system of Eqs. (16) is integrated with a 4th order Runge-Kutta algorithm with fixed time-stepping, and boundary conditions q = q = 0 to account for 10 0 0 N+1 the termination of the structure. The initial conditions -50 8 were set to zero. The time allowed for the elimination -100 of transients in each value of the frequency is 2000 T ( 0 n 6 -150 T = 2π/ω ), while the time interval τ over which the 0 0 0 averagepoweriscalculatedistypically1000T . Inorder -200 0 4 toconverttimeintervalstophysicalunitsonehastomul- -250 tiplybytheinverseoftheresonancefrequencyf−1,which 0 -300 2 gives τ 2.2 µs. In the linear regime, i.e., α = 0 and 0 ≃ -350 β =0,theaveragepowerP(dBm)densityplotonthesite number n - frequency f plane shown in Fig. 6 exhibits 0 0.5 1 1.5 2 -400 the expected behavior; that appreciable power transmis- f (GHz) sion only occurs for frequencies in the linear band, i.e., around f = 0.9 GHz. However, the activation of the nonlinear terms reveals that even a relatively low driv- FIG. 8: Average power density plot (in dBm) on the n−f ing amplitude ε = 0.3 is sufficient for nonlinear power plane, for λ = +0.01, γ = 0.01, α = −0.4, β = 0.08, and 0 transmission bands to bew formed. Thus, power can be ǫ0 =0.3. transmitted not only for frequencies in the linear band, butalsointheotherwiseforbiddenfrequencyregions. As When losses are small, the appearance of more non- it can be observed in Fig. 7, a number of transmission linear bands is possible, like the one seen just above the bands have appeared due to the nonlinearity. The ar- linear band in Fig. 7, at f 1 GHz. That transmission ∼ ray thus becomes transparent in power transmission at bandis due to aninterplaybetween nonlinearityandge- frequency intervals around the nonlinear resonances of a ometrical (lattice) resonances, and it is very sensitive to single varactor-loadedSRR oscillator. That self-induced changesintheparametersofthemodel. Notethatinthe transparency is a robust effect and survives with consid- density plots presented above the frequency is given in erably increased losses, as it is shown in Fig. 8, where natural units, using the relation f =f Ω=0.9Ω GHz. 0 8 It is very illuminating to see the time evolution of the 0 current flowing in each SRR of the array for some frequencies in the transmissionbands. We have chosentwo different frequencies, one in the linear band and one the ) m -50 high frequency band which is around the second har- B monic of the linear resonance frequency. The time evo- d lution of the currents in the first, fourth, seventh, and P( tenth (last) SRR are shown for f = 0.888 GHz and -100 f = 1.78 GHz in Figs. 9a and 9b, respectively. The currentsoscillateperiodicallyinallSRRs, althoughwith different amplitudes and relative phases with respect to 0 2 4 6 8 10 the driver. This figure suggests that power transmis- n sion for frequencies in the forbidden band gap occurs through the linear modes of the system. For frequen- FIG. 10: Power tranfer efficiency in different frequency re- cies outside the linear and nonlinear bands, the currents gions, i.e.,averagepower(indBm)asafunctionofsitenum- exhibit similar oscillations, althoughtheir amplitude be- ber n for an array ε0 = 0.3, λ = 0.01, γ = 0.001, N = 10, comes vanishingly small with increasing site number n. and f = 0.888 GHz (black diamonds); f = 1.197 GHz (red squares), f =1.78 GHz (green triangles). The relative efficiency of power transmission in the two the end of the extended array, and no reflection occurs. 2 The total number of SRRs is now N′ = N +N where 1 N =5,andthedissipationcoefficientγisnowafunction 1 n 0 i of the site number n, γ(n), given by -2 (a) γ , 1 n N γ(n)= γeµn , N≤<n≤ N′ , (34) (cid:26) ≤ 1 where µ = log(1/γ), with γ = 0.001 as before. For the N1 arraywithn dependentdissipationcoefficient,atypical − averagepowerdensityplotisshowninFig. 11,wherethe n 0 i power in each site is the average over a time interval of 5000 T . We observe, at least for the first N sites of the 0 (b) array,a pattern similar with those seen in Figs. 7 and 8. -1 That is, the appearance of the linear and the nonlinear 0 5 1τ0 15 20 transmission bands at about the same frequency intervals. We also observe that almost all of the transmitted power is absorbed from the five last sites of the array. FIG.9: Timedependenceofthecurrentin inthen−thSRR Importantly, the narrow band coming from geometrical for an array with N =10, λ=0.01, ǫ0 =0.3, γ =0.001, and resonancesinFigs. 7and8,whichisaroundf =1GHz, frequency (a) f = 0.888 GHz; (b) f = 1.78 GHz. The four different curves are for the current in the n=1 SRR (black- appearsin Fig. 11as well,and moreoverit seems that it solid), n = 4 SRR (red-dashed), n = 7 SRR (green-dotted), is capable of transmitting power to the end of the array. n=10 SRR(blue-dot-dashed). frequenciesusedinthe previousfigure,areshowninFig. 10, along with the power transmission for a specific fre- V. DISCUSSION AND CONCLUDING quency(f =1.197GHz)outsidealllinearandnonlinear REMARKS bands. It is remarkable that the transmission efficiency is almostthe sameforfrequencies inthe twopass-bands, We have used a simple electric circuit model of in- whilefortheotheronethetransmittedpowerpractically ductively coupled resistive-inductive-capacitive circuits vanish at the fourth site. with voltage-dependent capacitance to investigate the In the density plots above,the rightboundary wasac- transmission of power in a one-dimensional, end-driven, tually a reflectingone,allowingthe formationofstation- nonlinear magnetic metamaterial comprised of varactor- ary states in the array. We have also used a totally ab- loaded SRRs. We have used both a reflecting and an sorbing boundary by adding five more sites at the right absorbingboundaryfortheendthatisnotdrivenbythe end ofthe array,where the damping coefficientincreases external field. The importance of such discrete models exponentially with increasing site number. Thus, the fordescribingthedynamicbehaviourofSRR-basedmag- waves coming from the left are completely damped at netic metamaterials has been recently appreciated [27], 9 Their width apparently depends on the frequency, i.e., the bands at higher frequencies are significanlty wider 15 0 than those at lower frequencies. Moreover, we have also observed the formation, ei- -50 12 ther complete or partial, of another nonlinear pass-band -100 at frequencies slightly above the linear band (at f ∼ 9 -150 1 GHz). That band, which is rather narrow and can n -200 be observed both in Figs. 7 and 8 as well as in Fig. 11, 6 seem to depend sensitively on the model parameters λ, -250 f, ε (for α and β fixed), and mainly on the loss coef- 0 -300 3 ficient γ. Comparing Figs. 7 and 8, which have been -350 obtained for different loss coefficient γ =0.001 and 0.01, 0 0.5 1 1.5 2 -400 respectively, we see that the power transmission in that f (GHz) band considerably decreases when the losses increase by an order of magnitude. That type of nonlinear band appears due to geometrical resonances of the standing FIG.11: AveragepowerdensityplotPn(indBm)onthen−f waves that are formed in the array due to driving at plane for a varactor-loaded SRR array with absorbing right one end. Geometrical resonances seem to offer an im- boundary (see text), for λ = +0.01, ǫ0 = 0.297, γ = 0.001, portant mechanism for the generation of nonlinear pass- and N′ =15 (N =10). bands,whichcanbeemployedforconstructingtransmis- sion lines with banded power transmission spectra with bands at desired frequency ranges. Similar conduct-free althoughsimilardiscretemodelshavebeenemployedear- transmission lines could be formed by superconducting lierforthestudyoflocalizedexcitationsinthosesystems magnetic metamaterials, where the basic structural ele- [16, 35] ment(theSRR)isreplacedbyitsdirectsuperconducting In the present work, we have focused on the trans- analogue, the rf SQUID [17, 19]. The rf SQUID, where mission of power along the array, and especially that in the acronym stands for radio-frequency superconductive the nonlinear regime, which results in the formation of quantum interference device, is just a superconducting nonlinear pass-bands for frequencies close to the reso- ring interrupted by a Josephson junction, which makes nances of a single varactor-loadedSRR that is described it an intrinsically nonlinear element. A one-dimensional by the approximate model equation (13). The approx- arrayof rf SQUIDs coupled through their mutual induc- imate model holds fairly well for relatively low driving tances forms a nonlinear magnetoinductive transmission amplitudes,forwhichthedissipativediodecurrentofthe line[19]whosepowertransmissionpropertiesisamatter varactor can be neglected. It seems, thus that the prop- of future work. erties of the individual elements of the array predomi- nantly determine its transmission properties, at least in the case of weakmagnetoinductive coupling investigated Acknowledgments in the present work. The numerical results demonstrate that the nonlinear pass-bands transmit power with ef- ficiently comparable to that of the linear band. For the This workwassupportedin partby the EuropeanOffice parametersetsusedinthiswork,whichareclosetothose ofAerospaceResearchandDevelopment, AFSORaward in [50] for a single varactor-loaded SRR, there are three FA8655-10-1-3039, and by the EURYI and MEXT-CT- nonlinear bands, both above and below the linear band. 2006-039047. [1] Cui, W., Zhu, Y.-Y., Li, H.-x. & Liu, S. [2009] “Self- two-dimensional magnetic metamaterials,” Phys. Rev. E induced gap solitons in nonlinear magnetic metamateri- 80, 017601-1–017601-4. als,” Phys. Rev. E 80,036608-1–036608-5. [5] Enkrich,C.,[2005]“Magneticmetamaterials attelecom- [2] Cui, W., Zhu, Y.-Y., Li, H.-X. & Liu, S. [2010] “Soli- municationandvisiblefrequencies,”Phys. Rev.Lett.95, tonsexcitationsinaone-dimensionalnonlineardiatomic 203901-1–203901-4. chainofsplit-ringresonators,” Phys. Rev. E81,016604- [6] Feth, N., König, M., Husnik, M., Stannigel, K., Niege- 1–016604-9. mann, J., Busch, K., Wegener, M. & Linden, S. [2010] [3] Eleftheriou, M., Lazarides, N. & Tsironis, G. P. [2008] “Electromagneticinteractionofspit-ringresonators: The “Magnetoinductive breathers in metamaterials,” Phys. roleofseparationandrelativeorientation,”Opt. Express Rev. E 77, 036608-1–036608-13. 18, 6545–6554. [4] Eleftheriou,M.,Lazarides,N.,Tsironis,G.P.&Kivshar, [7] Geniet,F.&Leon,J.[2002]“Energytransmission inthe Yu. S. [2009] “Surface magnetoinductive breathers in forbiddenbandgapofanonlinearchain,”Phys.Rev.Lett 10 89, 134102-1–134102-4. rials: magnetism at optical frequencies,” IEEE J. Selec. [8] Gu, J., Singh, R., Tian, Z., Cao, W., Xing, Q., He, M.- Top. Quant. Electron. 12, 1097–1105. X.,Zhang,J.W., Han,J., Chen,H.&Zhang,W.[2010] [26] Litchinitser, N. M. & Shalaev, V. M. [2008] “Photonic “Terahertz superconductor metamaterial,” Appl. Phys. metamaterials,” Laser Phys. Lett. 5, 411–420. Lett. 97, 071102-1–071102-3. [27] Liu, N., Liu, H., Zhu,S. N. & Giessen, H.[2009] “Stere- [9] Hand,T.H.&Cummer,S.A.[2008]“Frequencytunable ometamaterials,” Nature Photon. 3, 157– . electromagnetic metamaterial using ferroelectric loaded [28] Maniadis, P., Kopidakis, G. & Aubry,S. [2006] “Energy split rings,” J. Appl. Phys. 103, 066105-1–066105-3. dissipation threshold and self-induced transparency in [10] Hesmer, F., Tatartschuk, E., Zhuromskyy, O., Rad- systems with discrete breathers,” Physica D 216, 121– kovskaya, A. A., Shamonin, M., Hao, T., Stevens, C. J., 135. Faulkner, G., Edwardds, D. J. & Shamonina, E. [2007] [29] Molina, M. I., Lazarides, N. & Tsironis, G. P. [2009] “Coupling mechanisms for split-ring resonators:Theory “Bulk and surface magnetoinductivebreathersin binary and experiment,” Phys. Stat. Sol. (B) 244, 1170–1175. metamaterials,” Phys. Rev. E 80, 046605-1–046605-9. [11] Huang,D.,Poutrina,E.&Smith,D.R.[2010]“Analysis [30] Pendry, J. B., Holden, A. J., Robbins, D. J. & Stewart, ofthepowerdependenttuningofavaractor-loadedmeta- W.J.[1999]“Magnetismfromconductorsandenchanced materialatmicrowavefrequencies,”Appl.Phys.Lett.96, nonlinear phenomena,” IEEE Trans. Microwave Theory 104104-1–104104-3. Tech. 47, 2075–2084. [12] Huang,Y.& McColl, W.F.[1997] “Analytical inversion [31] Penciu, R. S., Aydin, K., Kafesaki, M., Koschny, Th., ofgeneraltridiagonal matrices,” J. Phys. A: Math. Gen. Ozbay, E., Economou, E. N. & soukoulis, C. M. [2008] 30, 7919–7933. “Multi-gap individual and coupled split-ring resonator [13] Katsarakis, N., [2005] “Magnetic response of split-ring structures,” Opt. Express 16, 18131–18144. resonators in the far infrared frequency regime,” Opt. [32] Poutrina, E., Huang, D. & Smith, D. R. [2010] “Analy- Lett. 30, 1348–1350. sis of nonlinear electromagnetic metamaterials,” New J. [14] Kourakis,I.,Lazarides,N.&Tsironis,G.P.[2007]“Self- Phys. 12, 093010-1–093010-27. focusingandenvelopepulsegenerationinnonlinearmag- [33] Powell, D. A., Shadrivov, I. V., Kivshar, Yu. S. & neticmetamaterials,”Phys.Rev.E75,067601-1–067601- Gorkunov, M. V. [2007] “Self-tuning mechanisms of 4. nonlinear split-ring resonators,” Appl. Phys. Lett. 91, [15] Larouche, S., Rose, A., Poutrina, E., Huang, D. & 144107-1–144107-3. Smith, D. R. [2010] “Experimental determination of the [34] Shadrivov,I.V.,Morrison,S.K.&Kivshar,Yu.S.[2006] quadraticnonlinearmagneticsusceptibilityofavaractor- “Tunable split-ring resonators for nonlinear negative- loaded split ring resonator metamaterial,” Appl. Phys. index metamaterials,” Opt. Express 14, 9344–9349. Lett. 97, 011109-1–011109-3. [35] Shadrivov,I.V.,Zharov,A.A.,Zharova,N.A.&Kivshar [16] Lazarides, N., Eleftheriou, M. & Tsironis, G. P. [2006] Yu. S. [2006] “Nonlinear magnetoinductive waves and “Discrete breathers in nonlinear magnetic metamateri- domain walls in composite metamaterials,” Photonics als,” Phys. Rev. Lett. 97, 157406-1–157406-4. Nanostruct. Fundam. Appl. 4, 69–74. [17] Lazarides, N. & Tsironis, G. P. [2007] “rf superconduct- [36] Shadrivov, I. V., Reznik, A. N. & Kivshar Yu. S. ing quantum interference device metamaterials,” Appl. [2007] “Magnetoinductive waves in arrays of split-ring Phys. Lett. 90, 163501-1–163501-3. resonators,” Physica B 394, 180–183. [18] Lazarides, N., Tsironis, G. P. & Kivshar, Yu. S. [2008] [37] Shadrivov, I. V., Kozyrev, A. B., van der Weide, D. “Surface breathers in discrete magnetic metamaterials,” W. & Kivshar, Yu. S. [2008] “Tunable transmission and Phys. Rev. E 77, 065601(R)-1–065601(R)-4. harmonic generation in nonlinear metamaterials,” Appl. [19] Lazarides, N., Tsironis, G. P. & Eleftheriou, M. [2008] Phys. Lett. 93, 161903-1–161903-3. “DissipativediscretebreathersinrfSQUIDmetamateri- [38] Shalaev, V.M. [2007] , “Optical negative-index metama- als,” Nonlin. Phen. Compl. Syst. 11, 250–258. terials,” Nature Photonics 1, 41–48. [20] Lazarides, N., Molina, M. I. & Tsironis, G. P. [2009] [39] Shamonina,E.&Solymar,L.[2004]“Magneto-inductive “Breather induction by modulational instability in bi- wavessupportedbymetamaterialselements: components nary metamaterials,” Acta Phys. Pol. A 116, 635-637. foraone-dimensinalwaveguide,”J.Phys.D:Appl.Phys. [21] Lazarides, N. & Tsironis, G. P. [2010] “Driven linear 37, 362– . modes: Analytical solutions for finite discrete systems,” [40] Shelby, R. A., Smith, D. R. & Schultz, S. [2001] “Ex- Phys. Lett. A 374, 2179–2181. perimentalverificationofanegativeindexofrefraction,” [22] Lazarides, N., Molina, M. I. & Tsironis, G. P. [2010] Science 292, 77-79. “Breathersinone-dimensionalbinarymetamaterialmod- [41] Sievers,C.&Takeno,C.[1988]“Intrinsiclocalizedmodes els,” Physica B 405, 3007-3011. in anharmonic crystals,” Phys. Rev. Lett. 61, 970–973. [23] Lazarides, N., Molina, M. I., Tsironis, G. P. & Kivshar, [42] Smith,D.R.,Pendry,J.B.&Wiltshire,M.C.K.[2004] Yu. S. [2010] “Multistability and localization in coupled “Metamaterials and negative refractive index,” Science nonlinearsplit-ringresonators,”Phys.Lett.A374,2095- 305, 788–792. 2097. [43] Soukoulis,C.M.,Linden,S.&Wegener,M.[2007]“Neg- [24] Linden,S.,Enkrich,C.,Dolling,G.,Klein,M.W.,Zhou, ative refractive index at optical wavelengths,” Science J., Koschny, T., Soukoulis, C. M., Burger, S., Schmidt, 315, 47–49. F.&Wegener,M.[2004]“Magneticresponseofmetama- [44] Stevens,C.J.,Chan,C.W.T.Stamatis,K.&Edwards, terials at 100 Terahertz,” Science 306, 1351–1353. D.J.[2010]“Magneticmetamaterialsas1-Ddatatranfer [25] Linden, S., Enkrich, C., Dolling, G., Klein, M. W., channels: An application for magneto-inductive waves,” Zhou, J., Koschny, T., Soukoulis, C. M., Burger, S., IEEETrans. Microw.Theory Techniques58,1248–1256. Schmidt,F.&Wegener,M.[2006] “Photonicmetamate- [45] Sydoruk, O., Radkovskaya, A., Zhuromskyy, O., Sha-