Nonlinear magnetic metamaterials and possible applications on alloptical comparators and bistabilities in FabryPerot cavities
Abstract
We investigate the modulational instability and timedomain dynamics of nonlinear magnetic metamaterials composed of coupled splitring resonators loaded by Kerr nonlinearity. Our results indicate that the recently proposed optical switching of local optical index based on uniformresponse assumption seems fragile. We conceive two alternative schemes to utilize the valuable enhanced nonlinearity, one is to focus on fewbody systems and directly make use of the modulational instability (e.g., an optical comparator design), the other is to consider global switching arising from global feedbacks as in usual cases (e.g., FabryPerot cavities) rather than localresonancesbased switching of optical constants which may be destroyed by the discreteness and interactions as we show in this paper. We also try to provide comprehensive understanding of the relations between our results on discrete nonlinear metamterials and those on continuous metamaterials and natural materials in the literatures.
 PACS numbers

42.65.Pc, 42.65.An, 75.30.Kz, 64.60.F
I introduction
Linear metamaterials, ever since the theoretical works pendry1999 (); pendry1996 (); walser2001 (); pendry2000 () and the experimental works smith2000 (); shelby2001 (), have been highlighted by a great number of achievements. Nonlinear metamaterials, however, as Zheludev points out in [zheludev2010, ], are sill “green apples” on the knowledge tree for metamaterials.
In the discussions of nonlinear metamaterials, particularly those composed of nonlinear splitring resonators (SRR), it will be of great convenience to treat the metamaterials as continuous materials with controllable optical constants [e.g. and ], which allows us to employ existing concepts and methods in nonlinear optics. Several works since [zharov2003, ], including [lapine2003, ; agranovich2004, ; wen2006, ; joseph2010, ; roppo2010, ; zharova2005, ], focus on this approach. They either directly assume intensitydependentindex continuous media or assume uniform responses of the discrete system before calculating optical constants. These approaches share a common implicit assumption, i.e., on a scale much smaller than the wavelength, a cluster of discrete inclusions in a metamaterial should uniformly (at least slowlyvaryingly) respond to external radiation. This assumption justifies the continuous descriptions which investigate the materials by partial differential equations. Particularly, some of them ([zharov2003, ] and [chen2011, ]) propose applications for globalgeometryindependent optical switching based on the bistable properties of these local controllable optical constants. Note that these bistabilities zharov2003 (); chen2011 () are obtained after assuming that all the metaatoms respond uniformly in the materials, at least in the equal phase planes, just as in the linear cases.
However, we should note, that nonlinear metamaterials (usually composed of nonlinear SRR) are intrinsically discrete, and that they can even present strong local feedback at the “singleatom” level, unlike the conventional nonlinear optical materials to which the feedback is usually introduced globally (e.g. by two mirrors in a FabryPerot cavity). Therefore, we expect to gain more knowledge if we treat the nonlinear metamaterials as discrete, interacting nonlinear lattices as they are.
Some results have been obtained by treating the nonlinear metamaterials as discrete nonlinear lattices, including [lazarides2006, ; eleftheriou2008, ; cui2009, ; eleftheriou2009, ; lazarides2010, ; shadrivov2006, ; rosanov2012, ]. What we should pay particular attentions is the series of works lazarides2006 (); eleftheriou2008 (); eleftheriou2009 (); rosanov2012 (), in which the authors show that the steady responses of dissipative SRR arrays may violate the translational symmetry and can be in the multibreather states flach2008 () even when the arrays are driven by uniform radiation. The breather states in these works are theoretically constructed in the weekinteraction regime and in the anticontinuous limit aubry1997 () in which the states are a priori nonuniform and may in practice be excited by introducing chirping radiation lazarides2009 (); lazarides2009pre ().
In this paper, by investigating nonlinear effective circuit models, we show that even when we start from the nearlyzero states under uniform external radiations, modulational instabilities will step in the bistable ranges proposed in [zharov2003, ], i.e., the responses of the SRR array will spontaneously break the translational symmetry , and finally arrive at nonuniform states or even be chaotic, which will make the definition of optical constants in [zharov2003, ] invalid at least in the bistable range of interest. These results reveal that the switching of fielddependent optical constants may not exist due to modulational instabilities. Despite the negative results, we propose the possibility of constructing alloptical comparators based on the spontaneoussymmetrybreaking phenomenon of two strongly coupled SRR. We also propose to utilize the valuable enhanced nonlinearity in SRR arrays, by introducing global feedback mechanism as in usual cases. Finally, we comment on the relations of the results here based on our discrete models for metamaterials and those in the literatures.
Ii The theoretical model
When a plane wave is normally incident into a threedimensional nonlinear metamaterial composed of SRR, if we adopt the uniform assumption as in zharov2003 (), the SRR located on an equalphase plane of the wave will experience uniform external driving force (induced by other SRR and the external radiation) and respond with uniform amplitude and phase. We can thus simplify our considerations by focusing on a twodimension (or even onedimension) array of SRR radiated by uniform external fields. If we find that in some cases the uniform assumption is not valid even in the equal phase plane, we can undoubtedly conclude that such assumption is not valid in the threedimension bulk materials.
We model a rectangle lattice of SRR by LCR circuits coupled by mutual inductance (just as in zharov2003, ) with nearestneighbor approximations and periodic boundary conditions. In terms of dimensionless quantities, the dynamic equation of the lattice under slowlyvarying approximation is
where , are the (dimensionless) local electromotive force acting on the SRR at site and the electric charge hold by its capacitor, respectively. The is the dimensionless frequency; , the time; , the resistance; , mutual inductances between nearest neighbors in the ‘x, y’ directions, respectively. Then the steady states should satisfy the equation
(2) 
More detailed definitions of the quantities are presented in Appendix A.
Iii Modulational instabilities of uniform responses of an SRR array
The modulational instabilities of continuous nonlinear materials with negative refraction index have been studied in [wen2006, ], in which temporal, spatial, and temporalspatial modulational instabilities are all studied. The spatial modulational instabilities of threedimension discrete nonlinear SRR arrays induced by magnetoinductive waves have been studied in [shadrivov2006, ], in which it is shown that the unstable range may completely cover the decreasing range of the bistability, especially the two switching points. In this section, we in general employ the formalism in [shadrivov2006, ] but focus on the model (LABEL:dynamic) introduced in the previous section.
Under uniformresponse assumption as in [zharov2003, ], i.e., and , the Eqs.(2) for steady states are simplified
(3) 
We can solve this equation and obtain the bistabilities of uniform responses as illustrated in Fig.(1). However, since mutual interactions in densely stacked metamaterials are inevitably present, it is a nature question that whether possible excitations of magnetoinductive waves (spatial Fourier components of fluctuations) could induce modulational instabilities of these uniform states. By linear stability analysis (e.g.,[flach2008, ; shadrivov2006, ]), we can obtain the growth rate for each Fourier component of magnetoinductive waves Eq.(13). If for all Fourier components, has a negative real part, the corresponding state is stable. Otherwise, unstable.
After the linear stability analysis, the stable and instable uniform responses are denoted by solid and dashed segments in Fig.(1), respectively. Since we have chosen an equal phase plane to investigate, the mutual inductances between nearest neighbors should be positive in one direction and negative in the other direction. In all cases in Fig.(1), we can clearly identify critical points as the boundaries between the stable and instable ranges. As the system is tuned from the stable ranges (across the critical points) to the instable ranges, we expect that it will spontaneously break the translational symmetry and actually be in nonuniform states. We can also observe from Fig.(1) and Eq.(13) that in the low intensity limit, the modulational instability is absent for all cases implying that nonlinear metamaterials always allow uniform description in perturbative regimes.
Iv Timedomain response of a twodimension SRR array
A natural question arises that to what degree the actual response can deviate from the uniform assumption. In order to show an explicit picture of typical nonuniform responses when a system has evolved into the unstable range, we directly integrate Eqs.(LABEL:dynamic) using the RungeKutta method. In the recent papers lazarides2006 (); eleftheriou2008 (); eleftheriou2009 (), single breather states (localized states which break the translational symmetry) are constructed from the anticontinuous limit in the weakcoupling regime. But these methods possibly can not be used in our problems since the coupling may be too strong for a single breather’s existence. Here, we start from nearlyzero states with tiny random fluctuations and turn on the uniform radiation at the beginning. We then monitor the timedomain responses of twodimension model (containing 10 10 SRR) under different intensities of radiation. The results are shown in Fig.(2). The parameters are chosen as the case with uniform bistable response in Fig.(1). In that case, the critical radiation intensity is about which is the boundary between the stable range and unstable range. In Fig.(2), we test four radiation intensities , , and , of which the former three are near the critical value while the fourth is deep in the unstable range. The initial response is nearly zero with each individual response randomly distributing within . For each test intensity, we evolve the system for periods and plot the timedomain responses of the the SRR at two different sites, (1,1) and (5,5). The final states for each case are also plotted in 3D figures. We observe that for the intensity below the critical value , the response of the system quickly converge to a steady state which is uniform for all SRR, which means that the uniform assumption in [zharov2003, ] is justified. However, when the intensity enters the unstable range, e.g., , the system spontaneously evolve into nonuniform states and cannot even reach a steady state (it seems like a limit circle state with periodic oscillation). For the case very near to the boundary , the situation is similar to with the translational symmetry also drastically broken. When we come to the case in the deep unstable range , it is even worse is that the timedomain responses quickly become chaotic. A common characteristic for the cases with nonuniform responses is that the actual states are highly disordered with individual responses being randomly positive and negative, which can hardly support regular wave propagation.
V Optical comparators based on the spontaneous symmetry breaking
We can see from the above two sections that the responses of nonlinear metamaterials in the unstable ranges can be rather complex. If we want to make use of these nonuniform states, apart from directly investigating the details of the complexity (e.g., [cui2009, ; eleftheriou2008, ]) in these ranges, we may also focus on fewbody systems to avoid these complexities.
The spontaneous symmetry breaking (bifurcation) and its sensitivity on small deviations indicate that this phenomenon can be utilized to construct comparators. In Fig.3, we illustrate this possible application by plotting the evolution of the steady states (obtained from Eq.(2)) of two coupled SRR under the continuous tuning of the overall intensities of the two electromotive forces which deviate from each other by a small amount (). Although in the initial state the responses of the two SRR only have a tiny difference, the resulting states in the symmetrybreaking range can have a large contrast ratio. This property can be used to compare which signal (in our case the electromotive force) has a larger initial amplitude. Due to the small size of SRR compared to the working wavelength, the two SRR should be excited and detected via nearfield ways (e.g. [shadrivov2006oe, ]). Such a device, if the two SRR can be fabricated as identical as possible, can be used to realize ultrasensitive and possibly utlracompact alloptical comparators in applied areas. We should stress that such applications do not have to rely on SRR, and they may be even realized in any other devices if two identical nonlinear resonators are strongly coupled with each other, resulting possible spontaneous symmetry breaking (e.g., [davoyan2011, ; salgueiro2010, ; bulgakov2011, ]). These devices do not involve any assistance from electronic devices and thus may be of interest in the fields such as alloptical logic gates, alloptical signal processing and photonic computers (e.g., [mccutcheon2009, ; lu2011, ; wei2011, ]). One remarkable feature of our proposal is that we do not need high quality factor as other applications for optical signal processing mccutcheon2009 (); tanabe2005 (); kim2007 (); husko2009 (); hu2008 (), since we only assume in our calculations. We also do not need uniform bistability.
Vi Setting metamaterials into FabryPerot (FP) cavities
It is fascinating to see that due to field enhancement at the “singleatom” level, the nonlinearity of a single SRR can be strongly enhanced, and bistability can even occur for a single SRR at relatively low intensities, which certainly is an advantage of nonlinear metamaterials compared to natural materials. If SRR are simply stacked together, however, problems will arise as we have studied. One way out to utilize the enhanced nonlinearity is to focus on fewbody devices as we studied in the previous section, i.e, directly making use of the multistabilities and modulational instabilities. Another is to operate the stacked devices under the lowintensity radiation and introduce global feedback mechanism, e.g., setting them (e.g., the structure in [chen2011, ]) in FabryPerot cavities with optical lengths much larger than the wavelength. This approach can not only effectively employ the enhanced nonlinearity originating from the singatom level, but also avoid the modulational instabilities suffered by the devices at higher intensities. That is because the FPbistability requires a jump of the refraction index only of order (compared to the negativepositive jump proposed in [zharov2003, ; chen2011, ]), where is the wavelength and is the length of the cavity. Detailed estimations reveal that due to the field enhancement in the narrow gap, the required radiation intensity inside the FP cavity would be smaller by a factor of about than that required if the FP is filled with the same kind of pure Kerr materials (see Appendix C). This SRRenhanced nonlinearity may be even more powerful and useful when combined with other nonlinearityenhancing mechanisms, e.g., [hosono1992, ; shimizua1992, ].
Vii Discussions
We want to comment on the relations between nonlinear metamaterials composed of SRR and natural nonlinear materials composed of natural atoms. The usually mentioned bistabilities in natural materials are induced by a global feedback mechanism and work at a perturbative regime, i.e., the nonlinearity is restricted to the thirdorder. In other words, usual optical bistabilities are actually the nonlinear properties of certain global optical modes. No direct applications related to bistabilities in natural materials have been focusing on the nonperturbative regime in which the nonlinearity is very large at singleatom level or equivalently the radiation field intensity is comparable to the intrinsical field intensities in the atom level, since the responses may be very complex and very hard to control. Similarly, the proposals on controlling the bistabilities of the uniform optical constants zharov2003 (); chen2011 () are just examples of working in the nonperturbative regime of the effective index of the bulk metamaterials, so it is not at all surprising to see the unpleasant instability problems. Actually, the uniform bistabilities in [zharov2003, ; chen2011, ] are properties of local resonance responses rather global modes in usual case. Our proposal to introduce global feedback by setting an SRR array in an FP cavity is just an attempt to avoid these nonperturbative regimes and try to seek a perturbative regime (using optical modes) to utilize the SRRenhanced nonlinearity. We illustrate the difference between bistabilities of local responses and global modes in Fig.1.
We also want to stress on the relations between the spatial modulational instabilities discussed in this paper and those for continuous nonlinear media. Continuous descriptions assume homogenization of a material and employ partial differential equations. The modulational instabilities and related solitons in continuous descriptions (e.g., [wen2006, ]), if indeed valid, are actually triggered by certain exploding longwavelength (compared with the lattice constant) magnetoinductive waves shamonina2002 () (or equivalents in other discrete systems) and at the same time preserved by the decaying shortwavelength magnetoinductive waves in the discrete descriptions. However, our results in Fig.1 and Fig.2 show no indication of such longwavelengthonly modulational instabilities even in the low intensity range. All the modulationalinstabilityinduced nonuniform states in Fig.2 have shortwavelength components. This seems like a paradox. Remember that the selffocusinginduced spatial modulational instabilities originate from the propagations of waves. As we state in the theoretical model section, we have only analyzed the modulational instabilities of an equiphase surface and have not involved the stabilities of propagations of waves, and therefore our results do not imply the paradox. These arguments implies that at comparatively low intensities, the spatial modulational instabilities can still happen due to propagation effects. As an added note, the modulational instability discussed in this paper should be termed as subwavelength modulational instability as in a recent publication [noskov2012, ].
We should also note that in the perturbative range (or low intensity range) in Fig.1, the stable uniform responses enable the materials to be perturbatively characterized as nonlinear media with Kerr nonlinearity (also indicated in the Appendix C). Therefore, the works (mentioned in the Introduction section) which directly assume Kerr nonlinearities together with negative refraction indexes could be regarded as theories describing the perturbative properties of real discrete nonlinear metamaterials (e.g., as in lazarides2005 ()).
One may question our negative results by referring to the experimental works shadrivov2008 (); powell2009 (). Although frequency shift is demonstrated (certainly not sufficient for the index switching), no strong evidences (e.g., hysteresis curve) for bistabilityrelated index switching are discovered, however. Moreover, the nonlinear metamaterials in [shadrivov2008, ; powell2009, ] are excited by monopole antennas, a nonuniform experiment condition way too far from explicitly demonstrating the uniform bistability of optical constants. In fact, the two papers only claim that their metamaterials are tunable and do not claim that the bistability is confirmed. Another question would arise from the LCR model we employ. Although this model is more or less simplified and is not general, it can capture the main features of many metamaterials, e.g., resonances, loss, mutual interactions, Kerr nonlinearity. Therefore, our results can give strong caveat for considering the applications of the nonperturbative behaviors of effective optical constants (e.g., bistability of refraction index) when interactioninduced modulational instability in the lattice level is possibly present, since the effective optical constants in those cases may be ill defined. We also want to point out that the recent publication noskov2012 () mentioned above investigate another kind of nonlinear arrays. After making dipolar and slowlyvarying approximation near a resonance, they derive almost the same formulism as in this paper and also obtain modulational instability in the uniform bistable range.
We also feel it is necessary to clarify the difference between the bistability of optical constants in the proposals zharov2003 (); chen2011 () and the bistability induced by periodic structures herbert1993 (); xu1993 () in the literatures. The origin of the former is the local feedback at the singleatom level while that of the latter is just the periodicity, (termed as distributed feedback in the seminal paper winful1979 ()). The periodicity of the former has been homogenized, being replaced by a homogenized optical constant. What is more, the former is derived by local consistent equations of local optical constants and do not rely on the shape of the bulk material (explicitly shown in [zharov2003, ]; in [chen2011, ] different shapes are used with the same interpretation), while the latter is derived from the global consistent equations explicitly determined by the boundary or periodicity. One can also introduce such distributed or global feedbacks on a length scale larger than the wavelength to achieve bistability in metamaterials as we propose in the previous section. These differences all arise from the working wavelengths of these phenomena, one is larger than the lattice constant while the other is comparable or smaller than that. Therefore, the historical success herbert1993 (); xu1993 () cannot be evidences for the theoretical proposals in zharov2003 (); chen2011 ().
Finally, we should be reminded that if the whole array of SRR are rather sparsely arranged, i.e., the coupling between neighbors can be neglected (gas phase), the localbistabilitybased optical switching should also approximately be valid although the two switching points are inevitably affected according to Eq.(13). But in those case the array is so sparse that it can hardly be treated as effective materials. The coupling parameters chosen in this paper is more typical in realistic cases (condensed phase).
Viii Summary
By investigating the modulational instabilities and multistabilities, we find that some recently proposed optical switching based on the bistability arising from uniform assumptions may not be valid. In order to get rid of this dilemma, we propose two approaches to utilize the valuable enhanced nonlinearity of metamaterials at singleatom level, one is to focus on the multistabilities of fewbody problems (e.g., alloptical comparator) and the other is to operate them at low radiation intensities by introducing global feedback mechanism to achieve global bistability effects (e.g., FabryPerot cavity) in perturbative regimes of optical constants. Finally, we briefly discuss the relations between our results and existing ones in the literatures. Particularly, we distinguish the questioned bistability supported by local resonance and the usual bistability of global optical modes (global resonances, or global feedbacks) and stress that the former may easily be destroyed by modulational instability induced by discreteness and mutual interaction in dense materials.
Acknowledgement
We especially thank Cuicui Lu and Junhui Li for their careful reviews of this paper and constructive suggestions. We also thank the whole of Xiaoyong Hu’s group for helpful discussions and A. Alu for his comments and providing related literatures.
Appendix A The Theoretical model
The typical structure of SRR is shown in Pendry’s seminal paper pendry1999 (). According to the effective LCRcircuit model of the SRR in the static limit, we can write down the dynamic equation of SRR lattice.
(4) 
where , and is the effective selfinductance, resistance, and capacitance, respectively. The is the charge hold by the capacitor at the site . The are mutual inductances between nearest neighbors in the x (y) directions, respectively. Note that the nearestneighbor approximation is employed.
In our case where the split of the ring is filled by Kerr materials, the effective capacitance depends on the electric field intensity, , in the gap, i.e., . Here, the is the zerofield capacitance and the is the characteristic field intensity.
We also have the relation between and in the quasisteady states,
(5) 
where is the gap width of the capacitor. Then we arrive at the relation between and ,
(6) 
Here since we usually work in the weaknonlinearity regime of the Kerr material, the above relation can be approximated as,
(7) 
and
(8) 
Combining the Eqs. (4, 8) and make slowlyvaryingamplitude approximation, we can obtain the dynamic equation (LABEL:dynamic) after defining the following dimensionless quantities. , , , , , and .
For the convenience of the discussion in Sec.(VI), we write down the nonlinear magnetic susceptibility of a single SRR.
(9) 
where is the current and is the complex impedance. With the linear part separated, it reads
(10) 
Appendix B Linear stability analysis
We study the linear stability of uniform states solved from the steady state equation (3). We first introduce tiny fluctuations to the steady states, , and derive the linear dynamic equations for by use of the full dynamic equation.
(11) 
Since the lattice is uniformly periodic, the eigenmodes should be uniform plane waves. Thus let . We then have
where , and .The eigenvalues of the matrix on the left are the growth rates of the magnetoinductive fluctuations with wave vector .
(13) 
Appendix C Estimations in Sec.Vi
In the Sec.VI, we propose bistability by setting an SRR array into an FP cavity. First of all, for the bistability to occur the required radiation intensity should reach a certain scale to induce nonlinear deviation of the refraction index , where is the radiation wavelength in the metamaterial and is the length of the FP cavity. This is usually a small factor provided that , which will effectively reduced the required radiation intensity and thus help to avoid the modulational instabilities at higher radiation intensities. In this regime, the nonlinearity of the metamaterial can be treated perturbatively.
Thus we need that the radiation propagate at least distance before it significantly attenuate. To this end, we should work at the frequencies where the metamaterial is transparent, i.e., far from the absorption (, where the detuning ), and therefore the lossless approximation can be justified. According to the results in the Theoretical Model section, the effective magnetization coefficient of the metamaterials can be approximately written in this case as . To express in term of the radiation intensity , we should relate the dimensionless field intensity in the gap to the radiation intensity .
The total electromotive force accumulated along the ring is
(14) 
with being the radius of the SRR. Then in the lossless limit,the voltage across the gap will be
(15)  
Now, considering the relations and , we have by relating the above results,
Then , which indicates that selffocusing and selfdefocusing nonlinearity are both possible in the perturbative regime depending on the sign of the detuning and .
Now, we can estimate the radiation intensity required to induce the bistability of an FP cavity,
(for the bistability without FP cavity, the righthand side will be of .)then
Compared with an FP cavity filled with the same kind of pure Kerr material, the required intensity inside the cavity is reduced by a factor
(16) 
For the structure proposed in chen2011 (), if the can reach 100 (see Fig.3a of chen2011 ()), we can estimate . The expression for electromotive force (14) and the effective description should be still roughly valid for a wavelengthsize ratio as small as 5. The gapsize ratio is a tunable parameter depending on the fabrication technique which we assume can reach as small as 1/10. We also assume that the linear response (depending on the detuning and density of SRR) is not too small, i.e., . Then the reduction factor could be estimated to be about .
The role of the resonance can be assessed, from the above estimations, to be positive because it can render large field enhancement [ should better be small in Eq.(15)]; negative because it inevitably bring in energy loss that may make the material opaque ( should not be too small).
References
 (1) J. B. Pendry, A. J. Holden, D. J. Robbins, and W. J. Stewart, IEEE Transactions on Microwave Theory and Techniques 47, 2075 (1999).
 (2) J. B. Pendry, A. J. Holden, W. J. Stewart, and I. Youngs, Phy. Rev. Lett. 76, 4773 (1996).
 (3) R. M. Walser, A. P. Valanju, and P. M. Valanju, Phys. Rev. Lett. 87, 119701 (2001).
 (4) J. B. Pendry, Phy. Rev. Lett. 85, 3966 (2000).
 (5) D. R. Smith, W. J. Padilla, D. C. Vier, S. C. NematNasser, and S. Schultz, Phy. Rev. Lett. 84, 4184 (2000).
 (6) R. Shelby, D. Smith, and S. Schultz, Science 292, 77 (2001).
 (7) N. I. Zheludev, Science 328, 582 (2010).
 (8) A. A. Zharov, I. V. Shadrivov, and Y. S. Kivshar, Phy. Rev. Lett. 91, 037401 (2003).
 (9) M. Lapine, M. Gorkunov, and K. H. Ringhofer, Phy. Rev. E 67, 065601 (2003).
 (10) V. M. Agranovich, Y. R. Shen, R. H. Baughman, and A. A. Zakhidov, Phy. Rev. B 69, 165112 (2004).
 (11) S. Wen, Y. Wang, W. Su, Y. Xiang, X. Fu, and D. Fan, Phy. Rev. E 73, 036617 (2006).
 (12) A. Joseph and K. Porsezian, Phy. Rev. A 81, 023805 (2010).
 (13) V. Roppo, C. Ciracì, C. Cojocaru, and M. Scalora, J. Opt. Soc. Am. B 27, 1671 (2010).
 (14) N. A. Zharova, I. V. Shadrivov, A. A. Zharov, and Y. S. Kivshar, Opt. Expr. 13, 1291 (2005).
 (15) P.Y. Chen, M. Farhat, and A. Alù, Phys. Rev. Lett. 106, 105503 (2011).
 (16) N. Lazarides, M. Eleftheriou, and G. P. Tsironis, Phy. Rev. Lett. 97, 157406 (2006).
 (17) M. Eleftheriou, N. Lazarides, and G. P. Tsironis, Phy. Rev. E 77, 036608 (2008).
 (18) W. Cui, Y. Zhu, H. Li, and S. Liu, Phy. Rev. E 80, 036608 (2009).
 (19) M. Eleftheriou, N. Lazarides, G. P. Tsironis, and Y. S. Kivshar, Phy. Rev. E 80, 017601 (2009).
 (20) N. Lazarides, M. I. Molina, G. P. Tsironisa, and Y. S. Kivshar, Phy. Lett. A 374, 2095 (2010).
 (21) I. V. Shadrivov, A. A. Zharov, N. A. Zharova, and Y. S. Kivshar, Photonics and Nanostructures Fundamentals and Applications 4, 69 (2006).
 (22) N. N. Rosanov, N. V. Vysotina, A. N. Shatsev, A. S. Desyatnikov, and Y. S. Kivshar, Phys. Rev. Lett. 108, 133902 (2012).
 (23) S. Flach and A. V. Gorbach, Phys. Rep. 467, 1 (2008).
 (24) S. Aubry, Physica D 103, 201 (1997).
 (25) N. Lazarides, M. Molina, and G. Tsironis, Acta Physica Polonica A 116, 635 (2009).
 (26) M. I. Molina, N. Lazarides, and G. P. Tsironis, Phys. Rev. E 80, 046605 (2009).
 (27) I. V. Shadrivov, S. K. Morrison, and Y. S. Kivshar, Opt. Expr. 14, 9344 (2006).
 (28) A. R. Davoyan, I. V. Shadrivov, and Y. S. Kivshar, Opt. Lett. 36, 930 (2011).
 (29) J. R. Salgueiro and Y. S. Kivshar, Appl. Phys. Lett. 97, 081106 (2010).
 (30) E. Bulgakov and A. Sadreev, Phys. Rev. B 84, 155304 (2011).
 (31) M. W. McCutcheon, G. W. Rieger, J. F. Young, D. Dalacu, P. J. Poole, and R. L. Williams, Appl. Phys. Lett. 95, 221102 (2009).
 (32) C. Lu, X. Hu, Y. Zhang, Z. Li, X. Xu, H. Yang, and Q. Gong, Appl. Phys. Lett. 99, 051107 (2011).
 (33) H. Wei, Z. Wang, X. Tian, M. Kall, and H. Xu, Nature Communications 2, 387 (2011), doi: 10.1038/ncomms1388.
 (34) T. Tanabe, M. Notomi, S. Mitsugi, A. Shinya, and E. Kuramochi, Appl. Phys. Lett. 87, 151112 (2005).
 (35) M. K. Kim, I. K. Hwang, S. H. Kim, H. J. Chang, and Y. H. Lee, Appl. Phys. Lett. 90, 161118 (2007).
 (36) C. Husko, A. D. Rossi, S. Combrié, Q. V. Tran, F. Raineri, and C. W. Wong, Appl. Phys. Lett. 94, 021111 (2009).
 (37) X. Y. Hu, P. Jiang, C. Y. Ding, H. Yang, and Q. H. Gong, Nature Photonics 2, 185 (2008).
 (38) H. Hosono, Y. Abe, Y. L. Lee, T. Tokizaki, and A. Nakamura, Appl. Phys. Lett. 61, 2747 (1992).
 (39) A. Shimizu, M. KuwataGonokami, and H. Sakaki, Appl. Phys. Lett. 61, 399 (1992).
 (40) E. Shamonina, V. A. Kalinin, K. H. Ringhofer, and L. Solymara, J. Appl. Phys 92, 6252 (2002).
 (41) R. E. Noskov, P. A. Belov, and Y. S. Kivshar, Phys. Rev. Lett. 108, 093901 (2012).
 (42) N. Lazarides and G. P. Tsironis, Phys. Rev. E 71, 036614 (2005).
 (43) I. V. Shadrivov, A. B. Kozyrev, D. W. van der Weide, and Y. S. Kivshar, Opt. Expr. 16, 20266 (2008).
 (44) D. A. Powell, I. V. Shadrivov, and Y. S. Kivshar, Appl. Phys. Lett. 95, 084102 (2009).
 (45) C. J. Herbert and M. S. Malcuit, Opt. Lett. 18, 1783 (1993).
 (46) B. Xu and N.B. Ming, Phys. Rev. Lett. 71, 3959 (1993).
 (47) H. G. Winful, J. H. Marburger, and E. Garmire, Appl. Rev. Lett. 35, 379 (1979).