Theoretical description and design of nanomaterial slab waveguides: application to compensation of optical diffraction

Planar optical waveguides made of designable spatially dispersive nanomaterials can offer new capabilities for nanophotonic components. As an example, a thin slab waveguide can be designed to compensate for optical diffraction and provide divergence-free propagation for strongly focused optical beams. Optical signals in such waveguides can be transferred in narrow channels formed by the light itself. We introduce here a theoretical method for characterization and design of nanostructured waveguides taking into account their inherent spatial dispersion and anisotropy. Using the method, we design a diffraction-compensating slab waveguide that contains only a single layer of silver nanorods. The waveguide shows low propagation loss and broadband diffraction compensation, potentially allowing transfer of optical information at a THz rate. © 2018 Optical Society of America under the terms of the OSA Open Access Publishing Agreement OCIS codes: (160.4236) Nanomaterials; (260.2065) Effective medium theory; (230.7400) Waveguides, slab; (160.3918) Metamaterials; (310.6628) Subwavelength structures, nanostructures. References and links 1. P. Grahn, A. Shevchenko, andM. Kaivola, “Electric dipole-free interaction of visible light with pairs of subwavelengthsize silver particles,” Phys. Rev. B 86(3), 035419 (2012). 2. P. Grahn, A. Shevchenko, and M. Kaivola, “Electromagnetic multipole theory for optical nanomaterials,” New J. Phys. 14(9), 093033 (2012). 3. P. Grahn, A. Shevchenko, and M. Kaivola, “Interferometric description of optical metamaterials,” New J. Phys. 15(11), 113044 (2013). 4. P. Grahn, A. Shevchenko, and M., Kaivola, “Theoretical description of bifacial optical nanomaterials,” Opt. Express 21(20), 23471–23485 (2013). 5. A. Shevchenko, P. Grahn, and M. Kaivola, “Internally twisted spatially dispersive optical metamaterials,” J. Nanophoton. 8(1), 083074 (2014). 6. G. Shvets and I. Tsukerman, Plasmonics and Plasmonic Metamaterials, (World Scientific, 2012). 7. J. Ginn and I. Brener, “Realizing Optical Magnetism from Dielectric Metamaterial,” Phys. Rev. Lett. 108(9), 097402 (2012). 8. J. Valentine, S. Zhang, T. Zentgraf, E. Ulin-Avila, DA. Genov, G. Bartal, and X. Zhang, “Three-dimensional optical metamaterial with a negative refractive index,” Nature 455(7211), 376–379 (2008). 9. P. Moitra, Y. Yang, Z. Anderson, I. I. Kravchenko, D. P. Briggs, and J. Valentine, “Realization of an all-dielectric zero-index optical metamaterial,” Nature Photon. 7(10), 791–795 (2013). 10. S. B. Choi, D. J. Park, S. J. Byun, J. Kyoung, and S. W. Hwang, “Near-Zero Index: Optical Magnetic Mirror for Field Enhancement and Subwavelength Imaging Applications,” Adv. Opt. Mater. 3(12), 1719–1725 (2015). 11. J. B. Pendry, “Negative refraction makes a perfect lens,” Phys. Rev. Lett. 85(18), 3966 (2000). 12. I. I. Smolyaninov, Y. Hung, and C. C. Davis, “Magnifying Superlens in the Visible Frequency Range,” Science 315(5819), 1699–1701 (2007). 13. Z. Su, J. Yin, K. Song, Q. Lei, and X. Zhao, “Electrically controllable soft optical cloak based on gold nanorod fluids with epsilon-near-zero characteristic,” Opt. Express 24(6), 6021–6033 (2016). 14. N. Engheta, “Circuits with Light at Nanoscales: Optical Nanocircuits Inspired by Metamaterials,” Science 317(5845), 1698–1702 (2007). 15. D. A. B. Miller, “Optical interconnects to electronic chips,” Appl. Opt. 49(25), F59–F70 (2010). 16. D. A. B. Miller, “Are optical transistors the logical next step?” Nat. Photon. 4(1), 3–5 (2010). 17. M. Smit, J. Van Der Tol, and M. Hill, “Moore’s law in photonics,” Laser Photon. Rev. 6(1), 1-13 (2012). 18. A. F. Koenderink, A. Alu, and A. Polman, “Nanophotonics: Shrinking light-based technology,” Science 348(6234), 516–521 (2015). Vol. 26, No. 7 | 2 Apr 2018 | OPTICS EXPRESS 9134 #315543 https://doi.org/10.1364/OE.26.009134 Journal © 2018 Received 12 Dec 2017; revised 21 Mar 2018; accepted 22 Mar 2018; published 29 Mar 2018 19. D. K. Gramotnev and S. I. Bozhevolnyi, “Plasmonics beyond the diffraction limit,” Nat. Photon. 4(2), 83–91 (2010). 20. N. Kinsey, M. Ferrera, V. M. Shalaev, and A. Boltasseva, “Examining nanophotonics for integrated hybrid systems: a review of plasmonic interconnects and modulators using traditional and alternative materials,” J. Opt. Soc. Am. B 32(1), 121–142 (2015). 21. V. E. Babicheva, M. Y. Shalaginov, S. Ishii, A. Boltasseva, and A. V. Kildishev, “Long-range plasmonic waveguides with hyperbolic cladding,” Opt. Express 23(24), 31109–31119 (2015). 22. M. Silveirinha and N. Engheta, “Tunneling of Electromagnetic Energy through Subwavelength Channels and Bends using Epsilon-Near-Zero Materials,” Phys. Rev. Lett. 97(15), 157403 (2006). 23. Y. He, S. He, J. Gao, and X. Yang, “Nanoscale metamaterial optical waveguides with ultrahigh refractive indices,” J. Opt. Soc. Am. B. 29(9), 2559–2566 (2012). 24. T. Jiang, J. Zhao, and Y. Feng, “Stopping light by an air waveguide with anisotropic metamaterial cladding,” Opt. Express 17(1), 170–177 (2009). 25. Y. Huang, W. T. Lu, and S. Sridhar, “Nanowire waveguide made from extremely anisotropic metamaterials,” Phys. Rev. A 77(6), 063836 (2008). 26. T. Baba, “Slow light in photonic crystals,” Nat. Photon. 2(8), 465–473 (2008). 27. H. Kosaka, T. Kawashima, A. Tomita, M. Notomi, T. Tamamura, T. Sato, and S. Kawakami, “Self-collimating phenomena in photonic crystals,” Appl. Phys. Lett. 74(9), 1212–1214 (1999). 28. P. Belov and C. Simovski, “Canalization of subwavelength images by electromagnetic crystals,” Phys. Rev. B 71(19), 193105 (2005). 29. V. Kivijärvi, M. Nyman, A. Shevchenko, and M. Kaivola, “An optical metamaterial with simultaneously suppressed optical diffraction and surface reflection,” J.Opt. 18(3), 035103 (2016). 30. V. Kivijärvi, M. Nyman, A. Shevchenko, and M. Kaivola, “Optical-image transfer through a diffraction-compensating metamaterial,” Opt. Express 24(9), 9806–9815 (2016). 31. D. W. Prather, S. Shi, D. M. Pustai, C. Chen, S. Venkataraman, A. Sharkawy, G. J. Schneider, and J. Murakowski, “Dispersion-based optical routing in photonic crystals,” Opt. Lett. 29(1), 50–52 (2004). 32. D.W. Prather, S. Shi, J. Murakowski, G. J. Schneider, A. Sharkawy, C. Chen, B.Miao, and R.Martin, “Self-collimation in photonic crystal structures: a new paradigm for applications and device development,” J. Phys. D: Appl. Phys. 40(9), 2635 (2007). 33. D. R. Smith, D. C.Vier, T.H.Koschny, andC.M. Soukoulis, “Electromagnetic parameter retrieval from inhomogeneous metamaterials,” Phys. Rev. E 71(3), 036617 (2005). 34. C. Menzel, C. Rockstuhl, T. Paul, F. Lederer, and T. Pertsch, “Retrieving effective parameters for metamaterials at oblique incidence,” Phys. Rev. B 77(19), 195328 (2008). 35. D. Kwon, D. H. Werner, A. V. Kildishev, and V. M. Shalaev, “Material parameter retrieval procedure for general bi-isotropic metamaterials and its application to optical chiral negative-index metamaterial design,” Opt. Express 16(16), 11822–11829 (2008). 36. A. Andryieuski, R. Malureanu, and A. V. Lavrinenko, “Wave propagation retrieval method for chiral metamaterials,” Opt. Express 18(15), 15498–15503 (2010). 37. J. Yang, C. Sauvan, T. Paul, C. Rockstuhl, F. Lederer, and P. Lalanne, “Retrieving the effective parameters of metamaterials from the single interface scattering problem,” Appl. Phys. Lett. 97(6), 061102 (2010). 38. D. R. Smith, S. Schultz, P. Markos, and C. M. Soukoulis, “Determination of effective permittivity and permeability of metamaterials from reflection and transmission coefficients,” Phys. Rev. B 65(19), 195104 (2002). 39. A. Shevchenko, V. Kivijärvi, P. Grahn, M. Kaivola, and K. Lindfors, “Bifacial metasurface with quadrupole optical response,” Phys. Rev. Appl. 4(2), 024019 (2015). 40. A. Shevchenko, M. Nyman, V. Kivijärvi, and M. Kaivola, “Optical wave parameters for spatially dispersive and anisotropic nanomaterials,” Opt. Express 25(8), 8550–8562 (2017). 41. V. A. Markel, “Introduction to the Maxwell Garnett approximation: tutorial,” J. Opt. Soc. Am. A 33(7), 1244–1256 (2016). 42. A. Shevchenko, P. Grahn, V. Kivijärvi, M. Nyman, and M. Kaivola, “Spatially dispersive functional optical metamaterials,” J. Nanophoton. 9(1), 093097 (2015). 43. P. B. Johnson and R. W. Christy, “Optical constants of noble metals,” Phys. Rev. B. 6(12), 4370 (1972). 44. W. L. Bond, “Measurement of the Refractive Indices of Several Crystals,” J. Appl. Phys. 36(5), 1674–1677 (1965). 45. B. E. A. Saleh and M. C. Teich, , Fundamentals of Photonics (Wiley, 2007).


Introduction
Optical nanomaterials are often composed of specially designed nanoscale-size structural units. Nanomaterials that behave as homogeneous light-propagation media, are called optical metamaterials. In general, optical properties of these media are determined by local electromagnetic multipole excitations, and their optical response can therefore be tailored by shaping the structural units of the material [1][2][3][4][5]. In particular, higher-order multipoles, such as magnetic dipole and electric quadrupole can become prominent [6], giving rise to interesting optical phenomena. For example, optical metamaterials have been shown to exhibit macroscopic optical magnetism [7], counterintuitive negative refraction [8], and even more peculiar zero refractive index [9]. These unusual properties lead to novel applications, such as optical magnetic mirrors [10], perfect lenses [11,12], invisibility cloaks [13], and metamaterial-based components for Photonic Integrated Circuits (PICs) [14].
The development of PICs can help solve the energy and bandwidth problems of microelectronics [15][16][17][18]. However, in dielectric waveguides commonly used in photonics, the transverse confinement of light has a size scale determined by the diffraction limit, which makes PICs much larger than the corresponding electronic devices. This problem has led to the development of sub-diffraction limited plasmonic waveguides [19]. A typical drawback of these waveguides, however, is their high ohmic loss [20]. To simultaneously fulfill the criteria of sub-diffraction limited light confinement and low losses, some metamaterial waveguide designs have been proposed. They include slab waveguides with metamaterial claddings [21] and channels consisting of epsilon-near-zero metamaterials [22]. It is remarkable that indefinite-metamaterial waveguides with the refractive index reaching ultrahigh values can provide better confinement of light than metal waveguides [23]. In addition, metamaterial waveguides acting in a slow-light regime [24,25] can be used to enhance nonlinear light-matter interaction and create nonlinear photonic components, such as wavelength converters and optical buffers [26].
An interesting approach to light confinement and guidance consists of using materials that compensate for optical diffraction. In a diffraction-compensating medium, optical beams do not diverge upon propagation. Up to now, this effect has been demonstrated in photonic crystals [27], canalizing metamaterials made of effectively two-dimensional metal wires [28], and in metamaterials composed of three-dimensional nanoparticles [29,30]. Even though integrated circuits based on photonic crystals have been devised [31,32], the wavelength-scale size of the unit cells in these crystals limits the prospects of miniaturization. We have previously demonstrated that a metamaterial can be designed not only to suppress the diffraction effect, but also to significantly reduce the absorption and reflection losses exhibited by the structure [29], which facilitates integration of the material with other optical elements. To the best of our knowledge, however, diffraction-compensating metamaterial slab waveguides that can be suitable for use in PICs have not been studied even theoretically.
In this paper, we propose an approach to effectively design and characterize weakly absorbing nanomaterial slab waveguides that are inherently spatially dispersive and anisotropic. The structural units of the waveguides are assumed to be centrosymmetric but otherwise arbitrary. In analogy to three-dimensional metamaterials, nanomaterial slab waveguides can be characterized by the wave parameters, such as the refractive index and impedance, that depend on the wave propagation direction. However, in this case, the waves to be considered are the waveguide modes instead of plane waves. We call the new parameters the mode index n and impedance Z. Since the usual methods to retrieve n and Z [3-5, 33-39] cannot be applied to nanomaterial slab waveguides, we develop a new method that deals with waveguide modes and is based on their transmission and reflection coefficients. We use the method to design a nanomaterial waveguide with tailored spatial dispersion that provides diffraction-compensated propagation of the guided light. The waveguide contains only one layer of silver nanorods, which makes nanofabrication of the structure relatively straightforward. Moreover, the obtained diffraction-compensation effect is shown to be broadband, covering almost the entire fiber optical O band, within which the material exhibits low propagation losses.
The paper is divided into four sections. Section 2 presents the developed theory. Section 3 introduces the designed waveguide and describes its properties. Section 4 concludes the paper.

The method
The modes of a nanomaterial slab waveguide are infinite planar waves with the electric-field amplitude depending only on the coordinate along the normal to the slab [see the amplitude profile a(y) in Fig. 1(a)]. The mode profile is determined by the refractive indices n subs and n clad of the substrate and the cladding, respectively, and the structure of the nanomaterial (NM). The effective mode parameters n and Z that we introduce here are defined similarly to the wave parameters in bulk nanomaterials. The mode index is n = β/k 0 , where k 0 is the vacuum wave number and β is the mode propagation constant, and the impedance Z is determined by the ratio of the spatially averaged electric and magnetic fields of the mode. These parameters can be retrieved from the reflection and transmission coefficients of a segment of the waveguide [see the segment 2 in Fig. 1(a)] placed between two semi-infinite slab waveguides with a known core refractive index n 1 [segments 1 and 3 in Fig. 1(a)]. Any mode of segment 1 propagating through the nanomaterial strip will be partly transmitted and reflected. Thus, one can retrieve n and Z for this mode in the segment 2 from the reflection and transmission coefficients. The retrieval will be exact, if the mode profile a(y) is the same in the three segments.
To retrieve the mode parameters, we treat the segment 2 in Fig. 1(a) as a Fabry-Perot etalon shown in Fig. 1(b). This can be done, because at any fixed y, the waveguide modes behave approximately as plane waves, and consequently, their interface reflection and transmission coefficients can be written in terms of n and Z in analogy to Fresnel coefficients. The electric-field reflection and transmission coefficients, r and t, of the Fabry-Perot etalon are [40] Here, ρ i j and τ i j are the boundary reflection and transmission coefficients for a mode propagating from medium i to medium j, β z is the longitudinal propagation constant, and W is the width of the nanomaterial strip [see Fig. 1(b)]. The coefficients ρ i j and τ i j can be written as [40] where η k is a tangential impedance given by η k = Z k cos ϕ k for TE modes and η k = Z k /cos ϕ k for TM modes; ϕ k is the wave propagation angle with respect to z. The impedance Z 1 is the ratio of the transversely averaged electric and magnetic fields in segment 1, and Z 3 = Z 1 . The coefficients r and t are functions of the mode parameters n and Z, and therefore, these parameters can be solved from Eqs. (1) and (2). Numerical evaluation of r and t is explained later in this section.
For TE-polarized modes of segments 1 and 3, the mode index n in = β in /k 0 is connected to the corresponding core index n 1 by the equation Here, D is the thickness of the waveguide core [see Fig. 1(a)]. A similar equation is derived also for TM modes in Appendix A.
In general, the mode profile in segment 1 can differ from that in segment 2, which can lead to an additional optical loss that contributes to n and Z. However, often the retrieved mode parameters are very accurate even in the presence of the mode mismatch. Possible errors in the retrieval can be reduced if n 1 is iterated towards the value for which the mismatch is minimized. We use the following iterative method. First, n is retrieved from Eqs. (1) and (2). Next, the real part of the retrieved n is inserted into Eq. (5) as n in to obtain the corresponding waveguide core index n 1 . Then, the obtained core index replaces the previous n 1 and a new, more accurate value of n is retrieved. This procedure is repeated until the difference between two successively retrieved values of n is negligible. In the beginning, n 1 can be chosen intuitively, for example, equal to that of a bulk nanomaterial. If the structural units of the nanomaterial are much smaller than the wavelength, higher-order multipole excitations in them can be neglected and the effective medium parameters of the waveguide core can be calculated using the Maxwell-Garnett or Bruggeman approximations generalized to anisotropic materials [41].
The coefficients r and t are calculated numerically by considering the waveguide configuration of Fig. 1(a). In the calculations, we use the finite element method of the commercial software COMSOL Multiphysics. For each chosen value of n 1 , the mode profile is calculated analytically for the input waveguide (segment 1). Then, we generate this mode in COMSOL by introducing an electric-current distribution at a short distance from the front facet of the nanomaterial waveguide. The computation outputs are the field distributions in the segments 1 and 3. From these distributions, we obtain r and t in a way explained in the next paragraph. The mode generation by a distributed electric current is described in Appendix B.
In general, the incident mode has a wave vector β in =xβ x +ẑβ in,z [see Fig. 1(b)]; note that the tangential component β x is the same in all three segments. For TE modes, the electric field distributions related to the incident, reflected, and transmitted waves are of the form respectively. Hereû andv are unit vectors showing the direction of the electric field and b(y) is a unitless mode profile function taking values between −1 and 1 [similar to a(y) in Fig. 1(a)]. At normal incidence, we haveû =v. The quantities E i , E r , and E t represent the plane-wave like amplitudes of the incident, reflected, and transmitted modes, respectively, at the input and output facets of the nanomaterial waveguide. The coefficients r and t are As the computation directly yields the total electric field distribution E T (x, y, z) only, such that the coefficients r and t must be written in terms of E T (x, y, z). In principle, these coefficients can be calculated from single values of E T at a given coordinate y by using Eqs. (11) and (12). However, due to inevitable numerical errors, this can lead to some mistakes in the obtained values of r and t. In order to minimize the influence of these numerical errors, we average the results over the whole mode profile b(y) and obtain and Here, the fields are taken at distances d r and d t from the nanomaterial waveguide [see Fig. 1(b)]. Equations (13) and (14) are derived in Appendix C. These equations also yield the reflection and transmission coefficients for TM modes, if the electric-field related terms in the equations are replaced by the corresponding magnetic-field related terms. The validity of Eqs. (13) and (14) has been confirmed by direct numerical calculations. We note that r and t can be obtained by using any available computation methods, not only those based on COMSOL Multiphysics, or determined experimentally. Now that we have a method to calculate the mode parameters, we can describe any optical field propagating in the waveguide. We call such fields beams. Similarly to a plane-wave decomposition of a three-dimensional optical beam in free space, we can decompose a waveguide beam field into the waveguide modes by using the following inverse Fourier transform: where the mode amplitudeÊ(y, β x ) is found by Fourier transforming E(x, y, z). In the presence of optical anisotropy and/or spatial dispersion, the mode components can have different profile functions b(y, β x ). However, if the dependence of b(y, β x ) on β x is weak, which is often the case, we can factor the field functions as E(x, y, z) = b(y)E(x, z) andÊ(y, β x ) = b(y)Ê(β x ). Equation (15) then becomes Using Eqs. (15) and (16), we easily find the electric field of the beam at any point (x, y, z) as a superposition of the mode fields. The magnetic field can then be obtained by using the Maxwell-Faraday equation. For TM modes, a decomposition that is similar to Eq. (15) can be written for the magnetic field, from which the electric field can be calculated by using the Maxwell-Ampère equation.

A diffraction-compensating nanomaterial waveguide
Using the tools of Section 2, we have designed a diffraction-compensating nanomaterial slab waveguide, in which optical beams preserve their transverse intensity profiles I(x, y) upon propagation in the z-direction. In terms of waveguide modes, the effect requires that the longitudinal wave number β z does not depend on the mode propagation angle ϕ [see Fig. 1(b)] [29, 30, 42]. As β = k 0 n and β z = β cos ϕ, the diffraction-compensation condition is n(ϕ) = n(0)/cos ϕ. We note that n = n r + in i is a complex valued quantity. For waveguides of our interest, n i is small, and the condition is approximated to be n r (ϕ) = n r (0) cos ϕ .
To design a diffraction-compensating slab waveguide, we first retrieve n for angles ϕ = 0 and ϕ = ∆ϕ, where ∆ϕ can be chosen arbitrarily (for example, ∆ϕ = 15 • ). If the equality n r (∆ϕ) = n r (0)/cos ∆ϕ holds, we verify the diffraction-compensation condition also for other angles. If the equality does not hold, we adjust the nanomaterial structure. Furthermore, since n depends on the wavelength, the diffraction-compensation wavelength is found by calculating the whole spectra of n r (∆ϕ) and n r (0)/cos ∆ϕ and searching for a point at which the two curves cross each other, as explained later in connection with Fig. 3. The obtained wavelength can then be shifted by adjusting the design. The designed waveguide is shown in Fig. 2. It is composed of silver nanorods embedded in a 120 nm thick layer of ZnO on a SiO 2 substrate. The rods are 30 nm wide in the x and y directions and 130 nm long in the z direction. They form an array with periods of 80 nm in the x direction and 180 nm in the z direction. A similar structure was earlier used for a diffraction-compensating three-dimensional metamaterial [30]. To our surprise, the present structure allowed us to achieve diffraction compensation in a slab with only one layer of nanorods. This makes it possible to fabricate the waveguide using electron beam lithography and other standard nanofabrication methods. The diffraction-compensation in the material is based on the fact that the waveguide modes polarized in the xz plane couple to the longitudinal plasmon excitations of the nanorods (with the z-directed dipole moments) more efficiently than to the transverse ones (with the dipole moments perpendicular to z). The coupling is stronger at larger propagation angles ϕ with respect to the z-axis, leading to a higher mode index. The unit-cell dimensions are selected to   To verify that the diffraction-compensation condition is satisfied not only for ϕ = 15 • , we have calculated the functions n r (ϕ) and n i (ϕ) for the range of ϕ ∈ [−60 • , 60 • ]. Polar plots of these functions corresponding to λ = 1310 nm are shown in Fig. 4(a). Large angles are excluded, because the numerical errors in the applied finite element model become significant at these angles. The blue and red lines represent n r and tenfold magnified n i , respectively. The curve of n r is seen to be almost perfectly flat over the entire angular range between −60 • and 60 • . Therefore, the wavelength of 1310 nm is indeed the diffraction-compensation wavelength for beams that have their angular spectra within this range. The index n i is small at all propagation angles included in the plot, ensuring low absorption. Small peaks observed in the n i contour are associated with a slight mode mismatch between the nanomaterial waveguide and segment 1. The mismatch is caused by absorption in the nanorods and leads to an additional scattering that contributes to n i . At angles larger than ca. 20 • , the absorption is increased, but the profile mismatch due to the absorption is reduced. While the absorption is low, considering the plasmonic character of the diffraction compensation, we believe that it can be further reduced by the structure design. For example, the effect of diffraction compensation can possibly be achieved by using a high-index dielectric core and a structured cladding. The developed theory will still be applicable. An alternative way to calculate n r , which can also be used to verify the retrieved values of n r , is to calculate n r from the mode wavelength. To obtain the field periodicityλ z along z, the numerically calculated field distribution of the mode is fit with a decaying sinusoidal curve. The mode wavelength is thenλ =λ z cos ϕ and the index n r is calculated as λ/λ. Figure 4(b) shows, as an example, the fitting curve calculated for a mode propagating at an angle of 15 • . The fitted distribution of the x component of the field, E x (z), was taken at a height of ∆y = 120 nm from the waveguide surface, where it is already smooth. Inside the waveguide, the field distribution is heavily distorted by the nanostructures. The field E x (z) and its fit are represented by the black and green lines, respectively. The values of n r calculated in this way are shown as the white circles in Fig. 4(a). These values are in agreement with the results obtained from the reflection and transmission coefficients [the blue lines in Fig. 4(a)].
To obtain further confirmation of the diffraction compensation phenomenon in the designed waveguide, we used COMSOL Multiphysics to directly simulate the propagation of a strongly focused beam in a relatively large piece of the waveguide. Then we compared the simulation result with the beam propagation in a transparent homogeneous waveguide with n = n r (0). The beams were generated using surface-current distributions described in Appendix B with a Gaussian mode distribution in the x-direction. Figure 5(a) shows the longitudinal intensity profile I(x, z) of a Gaussian beam in the homogeneous and isotropic waveguide. The white dashed lines represent the radius of the beam defined as the distance from the beam axis, at which the intensity decreases by a factor of e 2 . Since the focal spot is extremely small, the beam strongly diverges. Its radius increases from 300 nm to about 1.5 µm within 1.6 µm of propagation from the beam waist. In the designed nanomaterial waveguide, the beam does not diverge at all within an even larger distance [see the centermost intensity profile in the Fig. 5(b)]. The two other beams shown in this figure have different wavelengths, 1280 nm and 1370 nm. Their divergence is also negligible. Therefore, the diffraction compensation phenomenon has a bandwidth of > 90 nm, corresponding to > 15 THz frequency, at which one can modulate the beam without affecting its diffraction-free propagation. For wavelengths well outside the considered range, the beam divergence does increase significantly.

Conclusions
We have developed an efficient theoretical approach that allows one to characterize nanomaterial slab waveguides in terms of the effective mode parameters, the refractive index and impedance. If the nanomaterial is spatially dispersive, the mode parameters depend on the propagation direction providing additional flexibility to the waveguide design. We have also developed a mode-decomposition method that makes it possible to evaluate electromagnetic field distributions of arbitrary beam-like fields propagating in the waveguide over arbitrarily large distances.
The mode parameters are retrieved from numerically calculated reflection and transmission coefficients of a short segment of the nanomaterial waveguide. In the calculations, the modes are generated by oscillating electric-current distributions. The reflection and transmission coefficients are derived from numerically calculated field distributions of the modes in a simple and convenient way. The introduced retrieval procedure is computationally light and easy to use for the design of nanomaterial slab waveguides.
Using the methods, we achieved a condition for diffraction-free propagation of optical beams in a metamaterial slab waveguide. The waveguide has a simple structure and can be fabricated by conventional nanofabrication methods. Surprisingly, only a single layer of metamolecules was found to be enough to provide diffraction compensation. The designed waveguide is not only conceptually new, but also promising for high-rate optical information transfer in photonic integrated circuits. The diffraction-compensation bandwidth covers the fiber optical O band, allowing tightly focused subpicosecond laser pulses to propagate in the waveguide without transverse spreading.

A. Derivation of the self-consistency equation
We assume that the field in the waveguide core is a superposition of two waves whose propagation paths are marked by the black and gray dashed lines in Fig. 6. For waveguides with n subs n clad , the intensity maximum of the wave superposition pattern is displaced from the middle of the core layer. To ensure that the vertical mode profile is preserved upon propagation, we write the following self-consistency equation for the two plane waves: where θ is the propagation angle of the waves, m is an integer, and φ 1 and φ 2 are the reflection phase shifts depicted in Fig. 6. The mode index n in is given by The phase shifts φ 1 and φ 2 can be found from where σ = 0 for TE polarization and σ = 1 for TM polarization. For TE polarization, the system of Eqs. (18)-(21) leads to the self-consistency condition given by Eq. (5). For TM polarization, the system leads to a condition given by tan k 0 n 2 1 − n 2 in D = − n 2 1 n 2 1 − n 2 in n 2 in − n 2 subs /n 2 subs + n 2 in − n 2 clad /n 2 clad n 2 1 − n 2 in − n 4 1 n 2 in − n 2 subs n 2 in − n 2 clad /(n 2 subs n 2 clad ) . (22)

B. Electric-current distribution as a source of a waveguide mode
The waveguide modes have the form E(x, y, z) = E 0 (y) exp(i(β x x + β z z)).
Consider an oscillating surface current in the z = 0 plane that radiates a mode of the above shape in the positive and negative z directions. Since the phase of the current distribution must match that of the mode, the current density can be written as The scalar product of this equation withv is r b(y) =v · E T (x, y, z) E i exp(−i[β x x − β in,z (z + W/2)]) −v ·û b(y) exp(i2β in,z (z + W/2)). (38) The field distribution E T is obtained by numerically solving the Maxwell equations for the given geometry and excitation mode. Since E T (x, y, z) includes near fields of the nanomaterial waveguide, we choose to calculate r from the field at a sufficiently large distance d r from the waveguide entrance facet. We then set z = −W/2 − d r in the above equations. We also choose x = 0. Multiplying Eq. (38) with b(y), integrating it over y ∈ [−∞, ∞], and dividing the result by the integral of b 2 (y), we obtain that is Eq. (13). Equation (14) is derived analogously. By integrating the fields E i , E r , and E t , we average them over the coordinate y, which significantly reduces the influence of numerical noise on the results. By multiplying the fields by b(y) before the integration, we calculate a weighted average, which reduces the dependence of r and t on the fields outside the waveguide, as well as on other possible modes, because their profiles are orthogonal to b(y). In addition, the multiplication of Eq. (38) with b(y) is necessary to have a non-zero integral for antisymmetric modes.