## Abstract

Microscopic structuring can change the effective properties of a material by several orders of magnitude. An example of this is animal bone, which has an effective elastic modulus that is more than 1,000 times larger than that of the constituent proteins. Here, we propose a broadband-enhancement principle of photonic nonlinearity that has a similar mathematical origin as the bone example. The proposed staggered array metamaterials violate the standard Miller’s rule in nonlinear optics and can enhance the third-order nonlinearity by more than a thousand to a billion times, depending on target operation frequencies. This metamaterial principle also enables manipulation of the individual components of the linear and nonlinear susceptibility tensors. Our biomimetic approach overcomes the fundamental speed-efficiency trade-off in current resonant enhancement schemes, making faster and more efficient all-optical devices possible for 1.55 μm wavelength. The principle is also applicable to ionic diffusion, heat conduction, or other transport problems.

## Introduction

Nonlinear photonic processes, such as harmonic generation, phase modulation, and multi-photon absorption, have diverse applications including laser sources and optical communications^{1,2,3}, three-dimensional lithography^{4}, bio-imaging^{5}, and many others^{6,7,8,9,10}. However, light–matter interaction is predominantly linear for common transparent photonic materials, thus necessitating the injection of high-intensity artificial light sources along the optimal crystallographic directions^{11}. Therefore, numerous studies have been conducted with the aim to attain photonic nonlinearity enhancement. Most of these studies have focused on the introduction of tailored electromagnetic resonances near the operation wavelength to boost light–matter interactions^{12,13,14,15,16,17,18,19}. Both quantum-electrodynamic approaches involving the engineering of intersubband transition energy levels^{12,13}, and classical electrodynamic enhancement methods—mostly utilizing photonic cavities or resonant structures^{13,14,15,16,17,18,19}—have been proposed. However, these resonant-enhancement schemes exhibit an inherent trade-off relationship between the degree of enhancement and the frequency bandwidth, thereby hindering their practical applications. Even for applications that do not require large bandwidths, large resonance-induced energy losses and a narrow fabrication error margin represent substantial obstacles. Moreover, only a few studies have investigated the controllability of the tensor form of nonlinear susceptibilities^{10,13,18}, despite the fact that it may have important practical consequences.

Here, we achieve both extreme enhancement and tensor form controllability of the linear and nonlinear susceptibilities of a nonlinear staggered array metamaterial (SAMM) over a broad bandwidth. The SAMM represents a periodic composite system of nonlinear materials and conductors with a deep subwavelength period. We analytically and numerically verify that the third-order nonlinear susceptibility can be enhanced by more than 10^{9} times at terahertz frequencies and by 10^{4} times at optical telecommunication frequencies. In addition, we show that 18 different crystal classes with second-order nonlinearity can be realized based on simple variations of the spatial configuration of a given nonlinear material and a conductor. The proposed broadband-enhancement principle is mathematically analogous to that of naturally occurring mechanical SAMMs, such as animal bones and bone-like materials (e.g., nacre and sea shell), which have enhanced elastic modulus by more than three orders of magnitude compared to their constituent soft proteins^{20,21}.

## Results

### Bespoke nonlinear susceptibility tensor of the proposed SAMM

The proposed SAMM structure is composed of alternating metal and nonlinear dielectric layers (Fig. 1a). The metal layers are composed of insulated metal patches in a planar lattice, and the nearest metal layers are staggered with half-unit-cell shifts in both the *x*- and *y*-directions (denoted as “A” and “B” layers in Fig. 1a). The metal and nonlinear dielectric layers (‘I’) are stacked in the sequence A–I–B–I–. The lateral unit-cell dimensions in the *x*- and *y*-directions, the lateral gap sizes in the same directions, and the thicknesses of the metal and dielectric layers are denoted as *a*_{x} and *a*_{y}, *g*_{x} and *g*_{y}, and *h*_{m} and *h*_{d}, respectively. Figure 1b depicts one unit cell of the SAMM.

The proposed microscopic (subwavelength) structure closely resembles the mechanical microscopic structure of mineralized collagen fibrils in animal bone and other hard biological materials (Fig. 1c)^{20,21}. In a simplified model, the fibrils in the bones are an array of isolated hard mineral plates embedded in a soft protein matrix, and the mineral plates are in a staggered formation as in Fig. 1a. Figure 1d schematically shows the microscopic deformation of the mechanical SAMM under macroscopic tensile strain. The mathematical parallelism between the enhancement of mechanical and electromagnetic properties is presented in the discussion section.

From the point of view of materials, SAMMs can be considered homogeneous, crystalline photonic materials because their periods are much smaller than the wavelength—as those of natural crystals—and their effective material properties, such as the linear and nonlinear susceptibilities and refractive indices, are well defined and remain constant regardless of the thickness of the sample^{22,23,24}. For these metamaterials, their linear (*n* = 1) and nonlinear (*n* > 1) homogenized electric susceptibilities \({\tilde{\mathbf{\chi }}}_{{\mathrm{eff}}}^{(n)}\) can be calculated based on the reciprocity theorem^{13,19,25} using the electric susceptibility tensors, \({\tilde{\mathbf{\chi }}}_{\mathrm{d}}^{(n)}({\mathbf{r}})\), of the constituent materials themselves and point-wise field-enhancement factors, \(F_{\alpha i}^{(\omega )}({\mathbf{r}}) = {\mathbf{\alpha }} \cdot {\mathbf{E}}_{{\mathrm{loc}}}^{(\omega )}({\mathbf{r}})/|{\mathbf{E}}_{{\mathrm{ave}}}^{(\omega )}|\) (with \({\mathbf{E}}_{{\mathrm{ave}}}^{(\omega )}||{\mathbf{i}}\)), which quantifies how strong the *α*-directional components of the local electric fields (\({\mathbf{E}}_{{\mathrm{loc}}}^{(\omega )}({\mathbf{r}})\)) are when the homogenized, macroscopic electric field (\({\mathbf{E}}_{{\mathrm{ave}}}^{(\omega )}\)) is along the *i*-direction. **α** and **i** are the unit vectors in the *α*- and *i*-directions, respectively (see Methods and Supplementary Note 1 for detailed derivation).

Using this approach, one can analytically derive the homogenized nonlinear susceptibility tensor of the proposed structure for second-harmonic generation (SHG). In this section, we assume for the sake of simplicity, the use of a metal with perfect electric conductivity and an insulator with frequency-independent nonlinear coefficients, thus resulting in the Kleinman symmetry condition^{11,26} formulated according to \(F_{\beta \gamma }^{(\omega )} = F_{\beta \gamma }^{(n\omega )} = F_{\beta \gamma }\). These restrictions are eliminated in the rigorous analysis of the actual enhancements and numerical demonstrations presented in the following sections. For the proposed structure, the local electric field is mainly oriented in the *z*-direction inside the dielectric layers owing to the parallel-plate capacitor-like geometry of the metal plates above and below the dielectric. Subsequently, the relevant field-enhancement factors are only *F*_{zx}, *F*_{zy}, and *F*_{zz}, and other factors are being negligibly small. The *z*-directional local electric field-enhancement factors under an *x*-directional macroscopic electric field, *F*_{zx}, of a specific structure (*a*_{x} = 4800 nm, *a*_{y} = 800 nm, *g*_{x,y} = 0.2*a*_{x,y}, and *h*_{d} = *h*_{m} = 20 nm) is plotted in Fig. 1e–f (see Methods). The actual profiles of *F*_{zx} and *F*_{zy} are very simple with a large and almost spatially uniform value (around 120 in the aforementioned case) throughout the space between the top and bottom metal plates, and with a negligible value elsewhere. These maximum field-enhancement values are approximately |*F*_{zx} | = 0.5*a*_{x}/*h*_{d} and |*F*_{zy} | = 0.5*a*_{y}/*h*_{d}^{27}. For the nonlinear dielectric, we assume that the microscopic pole directions are either the +*z*- or −*z*-direction or unpoled, thus resulting in the intrinsic second-order nonlinear susceptibility of the form \(\chi _{33}^{(2)}({\mathbf{r}}) = \chi _{33}^{(2)}p({\mathbf{r}})\) with *p*(**r**) = +1, −1, or 0, respectively. As a representative example, we consider a specific spatial pole configuration of *p*(**r**) = sign[*F*_{zy}(**r**)]. Figure 1g–h show *z*-directional local second-order nonlinear enhancement factor under an *x*-directional macroscopic electric field, \(\chi _{33}^{(2)}F_{zx}^2/|\chi _{33}^{(2)}|\), of this example. The nonlinear enhancement factor is very large (around 14,000) and spatially nearly uniform throughout the aforementioned high-field region between plates. The contracted notation^{11,26,28} of the effective nonlinear susceptibility simplifies to

if *g*/*a*_{x} and *g*/*a*_{y} are small (see Supplementary Note 2 for detailed derivation).

The tensor form in Eq. (1) corresponds to the monoclinic system with a single two-fold rotational symmetry. In a similar manner, one can realize various other nonlinear crystal systems by appropriately designing the spatial symmetry of *p*(**r**). If we also include the crystal designs in which the sides of the metal plates are intentionally arranged in a non-parallel orientation to any of the lattice vectors to break mirror symmetries, all of the noncentrosymmetric triclinic, monoclinic, orthorhombic, and tetragonal crystal classes can be realized (Fig. 2). If we consider the use of hexagonal plates in a hexagonal lattice, all of the noncentrosymmetric trigonal and hexagonal crystal classes are also covered (see Supplementary Fig. 1). Note that by tuning both the microscopic structure and the pole configuration in this simple manner, all noncentrosymmetric crystal classes—with the exception of the cubic classes—can be realized (a total of 18 point groups) and that only \(\chi _{zzz}^{(2)}\) (\(\chi _{33}^{(2)}\) in the contracted notation) of the constituent dielectric is utilized in all cases.

### Extreme enhancement of nonlinear susceptibilities at various frequency ranges

The local electric field-enhancement factors of the proposed SAMM show three noticeable characteristics. First, the magnitude of the enhancement factors can become very large for high-aspect-ratio structures (Fig. 1e–f, for example). The resulting nonlinear enhancements are even larger because they are products of multiple field-enhancement factors, depending on the order of the nonlinearity (Fig. 1g–h, compared to Fig. 1e–f). Additionally, these enhancement factors are nearly frequency-independent in the long-wavelength regime. A simple model, assuming perfect electric conductors and dispersion-less dielectrics, exhibits no frequency dispersion at all (*F*_{zx} = 0.5*a*_{x}/*h*_{d} and *F*_{zy} = 0.5*a*_{y}/*h*_{d}), and a rigorous model with realistic materials also exhibits a small frequency dependence, except for the cases with extremely high aspect ratios. Third, the enhanced local field appears over a significant portion (50% or even more) of the total volume of the effective medium, as opposed to small “hot spots”, as in previous tip-enhanced nonlinear structures. In structures with a hot spot near a sharp tip, or a narrow gap between metallic inclusions, the volume of the hot spots becomes smaller if one increases the degree of field enhancement by reducing the tip’s radius of curvature or the gap size. Subsequently, after homogenization with volume integration, the effective nonlinear susceptibility of these structures is not as drastically enhanced as their fields. These are the key factors that differentiate this non-resonant, polarization accumulation scheme of extreme nonlinearity enhancement via space-filling of local electric dipoles from many previous approaches. The conceptual explanation regarding to the polarization accumulation and the large susceptibility enhancement are in Supplementary Note 3 and Supplementary Fig. 2.

The effective third-order nonlinear electric susceptibility tensor, \({\tilde{\mathbf{\chi }}}_{{\mathrm{eff}}}^{(3)}\), has an overall enhancement factor which is proportional to the fourth power of *F*_{zx} or *F*_{zy}. The long-wavelength value of \({\tilde{\mathbf{\chi }}}_{{\mathrm{eff}}}^{(3)}\) for a SAMM with *a*_{x} = *a*_{y} = *a* becomes

based on the simple model, and it satisfies the tetragonal symmetry (4*mm*) (see Supplementary Note 2 for detailed derivation). In the following section, we focus on \(\chi _{{\mathrm{eff}},11}^{(3)}\), which leads to an extremely large third-order nonlinear polarization along the *x*-direction under an incident *x*-directional electric field.

For more accurate calculations in the rest of the paper, we now adopt a rigorous analytical model, considering the actual complex permittivity of both the metal and the dielectric. This model is used to calculate the effective refractive index and the nonlinear susceptibilities (see Supplementary Note 4 for detailed derivation). Fig. 3 shows the analytically calculated and numerically simulated \(\chi _{{\mathrm{eff}},11}^{(3)}\), over several decades of frequencies (see Methods). In Fig. 3a, the calculated enhancement factor \(\chi _{{\mathrm{eff}},11}^{(3)}/\chi _{33}^{(3)}\) from the rigorous model is depicted as a function of both the aspect ratio of the SAMM, *a*/*h*_{d}, and the frequency, with the structural parameters of *h*_{d} = *h*_{m} = 10 nm and *g* = 0.2*a*. The lateral unit-cell dimension, *a*, is changed to vary the aspect ratio. An almost uniform, extreme enhancement is achieved from near-zero frequency up to the lowest resonance frequency. For low frequencies, the simple estimation of \(\chi _{{\mathrm{eff}},11}^{(3)}\) in Eq. (2) matches the rigorous model well. As the frequency increases, material-induced resonances occur when \(a^2/(4h_{\mathrm{m}}h_{\mathrm{d}})\sim |\varepsilon _{\mathrm{m}}^{(3\omega )}/\varepsilon _{\mathrm{d}}^{(3\omega )}|\) (gray dashed line) and when \(a^2/(4h_{\mathrm{m}}h_{\mathrm{d}})\sim |\varepsilon _{\mathrm{m}}^{(1\omega )}/\varepsilon _{\mathrm{d}}^{(1\omega )}|\) (black dashed line) (see Supplementary Note 4 for details about the resonance conditions). These resonances can be classified as electrostatic resonances, with a resonance condition determined by quasi-static boundary conditions^{29,30,31,32,33}. Because the magnitude of the metal permittivity is reduced as the frequency increases, the resonance frequency can be increased by lowering the aspect ratio. Therefore, one should design the structural parameters by considering the target frequency range, the amount of enhancement desired, and the acceptable amount of dispersion. To reduce the computational load, we consider structures that are uniform in the *y*-direction (*g*_{y} = 0) in the following calculations. \(\chi _{{\mathrm{eff}},11}^{(3)}\) for non-zero *g*_{y} is also of the same order of magnitude as long as *g*_{y}/*a*_{y} ≪ 1.

Figure 3b shows the enhancement factor (\(|\chi _{{\mathrm{eff}},11}^{(3)}/\chi _{33}^{(3)}|\)) for two different aspect ratios calculated by a nonlinear transfer-matrix-based parameter extraction method using nonlinear finite-difference time-domain (FDTD) simulation results. The values are compared to the analytical model predictions. For a moderate aspect ratio (*a* = 900 nm, *g* = *a*/5, *h*_{d} = *h*_{m} = 10 nm, i.e. *a*/*h*_{d} = 90), we observed a million-fold enhancement of the effective nonlinear susceptibility compared to that of the constituent nonlinear dielectric (Fig. 3b, blue line), which is already two orders of magnitude higher than that of a previously proposed theoretical resonant structure^{17}. Moreover, the enhancement is almost uniform from a near-zero frequency up to 2 THz, which is much lower than the electrostatic resonance frequency (5 THz). A SAMM with a larger aspect ratio (*a* = 6000 nm, *g* = *a*/5, *h*_{d} = *h*_{m} = 10 nm, i.e. *a*/*h*_{d} = 600) shows more than a billion-fold enhancement (Fig. 3b, green line). Although there is an electrostatic resonance at a sub-terahertz frequency for the structure with this high aspect ratio, the enhancement remains extremely large from the zero to the terahertz frequency range. The enhancement factors from the FDTD simulations and the rigorous analytical model are in excellent agreement. The retrieved \(\chi _{{\mathrm{eff}},11}^{(3)}\) value remained constant regardless of the thickness and the input intensity (see Supplementary Fig. 3). Therefore, the extracted \(\chi _{{\mathrm{eff}},11}^{(3)}\) is well defined as an effective bulk material property of the SAMMs.

At higher frequencies, including the infrared and visible, the highest aspect ratio that one can use for broadband enhancement becomes lower owing to the smaller magnitude of the metal permittivities. Figure 3c shows the enhanced third-order nonlinearity of the proposed structure with *a* = 75 nm, *g* = 15 nm, and *h*_{d} = *h*_{m} = 5 nm with realistic materials at 200 THz, which is close to an important optical-fiber communication band. The electrostatic resonances occur at 67 THz and 200 THz. The \(\chi _{{\mathrm{eff}},11}^{(3)}\) enhancement is larger than 4×10^{4} at 200 THz. In addition, the enhancement in the long-wavelength limit is larger than 10^{3}, which is still a very high enhancement.

### Thin and efficient nonlinear reflectors utilizing SAMMs incorporating nonlinear quantum wells

The proposed microscopic structure of the staggered conductors functions as a universal magnifier and symmetry-modifier of the intrinsic nonlinear susceptibility tensor of the constituent nonlinear material of choice, regardless of the latter. Accordingly, the quantitative enhancement factors shown in Eqs. (1) and (2) as well as the symmetry-changing capabilities demonstrated in Fig. 2 remain the same when we change the nonlinear materials. Thus, one can utilize quantum wells (whose sub-band energy levels are carefully designed to possess resonant-enhanced nonlinear susceptibilities that are several orders of magnitude larger than those of typical nonlinear materials in their off-resonance regime) as the nonlinear material between the conductors and the microscopic structure will further enhance these large nonlinearities, by another several orders of magnitude. This allows the use of metasurface types of nonlinear reflectors for massively parallel, pixel-by-pixel nonlinear optical operations (Fig. 4), rather than using waveguides or crystal-type devices with long propagation paths.

The nonlinear reflectors in Fig. 4a comprise a thin SAMM slab with GaN/AlN quantum well^{34}, dielectric impedance-matching layers above and below^{35,36}, and a metallic mirror (see Methods). The electromagnetic properties of the GaN/AlN quantum well is summarized in Supplementary Fig. 4.

Figure 3b shows the switching performance of the optimized nonlinear reflector (*a* = 57 nm, *g* = 20 nm, *h*_{d} = 5.1 nm, *h*_{m}/2 = 10.2 nm, *h*_{u} = 357 nm, and *h*_{l} = 57 nm) for pulses which are 30 fs with a center wavelength of 1.55 μm. The thickness of the SAMM and the total thickness including dielectric layers are only of the order of 1/60 and 1/3 of the vacuum wavelength, respectively. While resonant-cavity-based all-optical switching approaches sacrifice speed and bandwidth in exchange for lower switching power[2], our non-resonant-enhancement scheme does not suffer from this limitation, which is apparent from the 26 fs output pulses (red curve in Fig. 4b inset). The recovery speed is ultimately limited by the relaxation time of the intersubband transitions in quantum wells, which can be of the order of 100 fs or faster^{34,37}. The system exhibits a modulation in excess of 10 dB with a pulse energy utilization of only 380 fJ, while the reflectance is more than 50% at the on-state mode of operation. This 10 dB modulation pulse energy is a record-low value for all-optical switches with similar modulation speed[2], and is 250 times smaller than previous quantum-well-based all-optical switches^{37}. We visually compare switching performance of the result in Fig. 4 with several state-of-the-art all-optical switches^{37,38,39,40,41,42,43} in Supplementary Fig. 5 using a Pareto optimality plot based on the switching time and the required switching energy. The maximum local electric field magnitudes are of the order of 5 × 10^{8} V m^{−1} for the 30 fs, 380 fJ pulses, which are below the damage threshold for many solid dielectric materials^{44}. Although the effective susceptibilities are now dispersive owing to the resonant nature of the quantum wells, the bandwidth is still large enough to cover a significant part of the optical communication bands or to allow femtosecond-pulse operations. In addition, we verify that the designed device maintains more than 5 dB modulation depth with same pulse energy (380 fJ) under up to 30% (±10 nm) lateral miss-alignment (zero miss-alignment being the perfectly staggered case) between two metal layers of the SAMM slab, relative to the lateral unit cell size (Fig. 4c).

## Discussion

Nonlinear susceptibilities of many natural nonlinear crystals can be estimated using Miller’s rule, which states that the *n*th-order nonlinear susceptibilities (e.g., \(\chi ^{(2)}(2\omega ;\omega ,\omega )\)) can be approximated to an (*n* + 1)th-order linear susceptibilities (e.g., \(\chi ^{(1)}(2\omega )\chi ^{(1)}(\omega )\chi ^{(1)}(\omega )\)) with a constant factor, *δ*, known as the Miller’s delta^{11,19,45}. However, according to O’Brien et al.^{19}, the nonlinear response predicted by Miller’s rule can be several times larger or smaller than the measured values for composite media such as metamaterials. In case of the proposed SAMMs, the actual nonlinear susceptibilities (Eqs. (1) and (2)) can differ from those predicted by Miller’s rule, by several orders of magnitude. Furthermore, the factor of discrepancy can be expressed as [*M*_{d}]^{−n}[*M*_{s}]^{−1}, where *M*_{d} and *M*_{s} are the enhancement factors of the effective oscillator density and strength, respectively, in a homogenized oscillator model of SAMMs (see Supplementary Note 5 for detailed derivation). Given that *M*_{d} and *M*_{s} are *a*/(4*h*_{d}) and *a*/(2*h*_{d}), respectively, for the proposed SAMMs^{27} and can be larger than 1000, the predictions using the original Miller’s rule without these correction factors may be 10^{9} and 10^{12} times different from the true values for \(\chi _{{\mathrm{eff}}}^{(2)}\) and \(\chi _{{\mathrm{eff}}}^{(3)}\), respectively. In Fig. 5, we compare the 3rd order nonlinear susceptibility enhancement factors (\(\chi _{{\mathrm{eff}}}^{(3)}/\chi _{\mathrm{d}}^{(3)}\)) of the SAMMs in this paper (blue crosses) with the estimated values based on the conventional Miller’s rule (black dots) and the modified Miller’s rule developed in this study (red dots) (see Methods). It is confirmed that the modified Miller’s rule with \(\delta^{\prime} = \delta [M_{\mathrm{d}}]^{ - 3}[M_{\mathrm{s}}]^{ - 1}\) estimates third-order nonlinear susceptibilities of SAMMs well. While large effective oscillator strength enhancement factor *M*_{s} can be found in many previous studies^{16,17,31,33,46,47}, large effective oscillator density enhancement factor *M*_{d} is a unique characteristic to SAMMs and induces a dramatic deviation from the conventional Miller’s rule. In fact, our SAMMs can be considered as an example of a new category of electromagnetic media with designable oscillator densities (see Supplementary Note 5 for details).

The enhancement of the electric susceptibilities demonstrated herein can be viewed as a special case of a general enhancement principle of effective material properties that occur in two-phase composite materials. For example, the enhancement of the elastic modulus of the mineralized collagen fibrils in bone (Fig. 1c, d) can be conceptually explained as follows: First, the macroscopic tensile strain induces a shear strain enhancement by a factor of 0.5*a*/*h*_{d} inside the protein layer, which results in an enhanced shear stress by a factor of 0.5*a*/*h*_{d}. This enhanced shear stress accumulates over the distance of 0.5*a* and appears as tensional stress in the mineral plates with a thickness of *h*_{m}, thereby inducing another enhancement factor equal to 0.5*a*/*h*_{m}. Rigorous derivations of the effective mechanical properties can be found in previous reports^{20,21}, and the expression of the effective elastic modulus is similar to Eq. (3) in Methods with *n* = 1 (linear case) and becomes *M*_{eff} ≈ *G*_{p}(*a*/*h*_{d})^{2}/8, where *G*_{p} is the shear modulus of the protein matrix. This shows that the system can have an effective elastic modulus that is several orders of magnitude larger than that of the host material, for cases involving high aspect ratios (*a*/*h*_{d} ≫ 1). This mechanical enhancement principle has been previously adopted for biomimetic hard and tough materials made of ceramics and polymers^{48}. We note that this mechanism for enhanced elastic modulus is analogous to the explanation in Supplementary Note 3.

In general, we argue that the same principle can be encountered in various homogenization problems (including the examples already discussed) with a constituent equation of the form **F** = *σ***E**, where **F** is the relevant flux density (e.g., heat or mass), **E** is the conservative driving force (e.g., the temperature gradient or mass density gradient), and *σ* is the conductivity. **E** is an irrotational vector field owing to the conservative nature, and **F** becomes divergence-free in a steady state. Subject to these conditions, the effective conductivity for a two-phase composite material with isolated high-conductivity phases embedded in a continuous low-conductivity phase host material in a staggered formation, as in our SAMM, is *σ*_{eff} = *σ*(*a*/*h*_{d})^{2}/8 (see Supplementary Note 6 and Supplementary Fig. 6 for detailed derivation). Thus, the effective properties of the composite materials can be many orders of magnitude higher than those of the continuous host, which may occupy more than half of the total volume, while the enhancement factor is purely geometrically controlled. Moreover, the nonlinearity of the conductivity is enhanced by even higher orders of magnitude in a manner identical to that of the enhancement of the nonlinear susceptibility, which might find use in nonlinear control of heat and mass flow, for example.

It has been shown that the symmetry and strength of individual components of linear and nonlinear susceptibility tensors can be tuned by the appropriate design of microscopic composites of nonlinear dielectrics and metals and that 18 nonlinear crystal classes can be realized. A near dispersion-less enhancement of *χ*^{(3)} of the order of 10^{6} was demonstrated for broad terahertz frequency ranges, and resonant enhancements of the orders of 10^{9} and 10^{4} were demonstrated for terahertz and near-infrared frequencies, respectively. Ultrafast, low-power, pixel-by-pixel all-optical modulation schemes at optical communication wavelengths was presented using a structure with a subwavelength thickness (*λ*/3). These SAMMs may serve as a useful platform for scientific studies of photonic nonlinearity with weak light or for high-density integration of nonlinear photonic devices. Moreover, the enhancement principle may find applications in other branches of science, while synergistic enhancements of multi-physical interactions are also imaginable.

## Methods

### Linear and nonlinear effective electric susceptibility tensors of metamaterials

The general form of the effective electric susceptibility tensors for a three-dimensional metamaterial can be formulated in a similar manner to that of the two-dimensional case (metasurface), which was previously derived from the Lorentz reciprocity theorem^{13,19,25}. In general, the *n*th-order nonlinear susceptibilities are (*n* + 1)th-order tensors. However, for *n*th-harmonic generation problems, the corresponding nonlinear susceptibilities can be represented with a more concise second-order tensor form using a contracted notation^{11,26,28},

where *V* is the unit-cell volume and \({\tilde{\mathbf{{\chi }}}}_{{\mathrm{eff}}}^{(n)}\) is a \(3 \times \left( {\left( \begin{array}{c}3\\ n\end{array} \right)} \right)\) homogenized nonlinear susceptibility tensor where \(\left( {\left( \begin{array}{c}3\\ n\end{array} \right)} \right) = (2 + n)!n!^{ - 1}2^{ - 1}\) is a combination with repetitions. \({\tilde{\mathbf{F}}}_n^{(\omega )}({\mathbf{r}})\) is a \(\left( {\left( \begin{array}{c}3\\ n\end{array} \right)} \right)\)-by-\(\left( {\left( \begin{array}{c}3\\ n\end{array} \right)} \right)\) tensor of field enhancements. Each component of \({\tilde{\mathbf{F}}}_n^{(\omega )}({\mathbf{r}})\) is a combinatorial summation of the products of *n* individual field-enhancement factors, \({E}_{{\alpha}i}^{({{\omega})}}({\bf{r}})={\bf{\alpha}}\cdot{\bf{E}}_{{\mathrm{loc}}}^{({{\omega})}}({\bf{r}}) / |{\bf{E}}_{{\mathrm{ave}}}^{({{\omega})}}|\) (with \({\mathbf{E}}_{{\mathrm{ave}}}^{(\omega)}{\parallel} \, {\mathbf{i}}\)), which quantifies how strong the *α*-directional components of the local electric fields (\({\mathbf{E}}_{{\mathrm{loc}}}^{(\omega )}({\mathbf{r}})\)) are when the macroscopic electric field (\({\mathbf{E}}_{{\mathrm{ave}}}^{(\omega )}\)) is in the *i*-direction. \(\tilde{\chi}^{(n)}{(\bf{r})}\) is the local \(3 \times \left( {\left( \begin{array}{c}3\\ n\end{array} \right)} \right)\) nonlinear susceptibility tensor at location **r** for the *n*th harmonic generation (for *n* = 1, Eq. (3) is reduced to the homogenized linear susceptibility). The derivation and explicit, element-wise form of Eq. (3) is given in Supplementary Note 1.

### Local electric-field-enhancement factor plot in Fig. 1

A lossy Drude model is assumed for gold, which is the metallic material (plasma frequency: 2.18 PHz, damping frequency: 6.45 THz)^{49}. Crystalline zinc oxide is assumed to be the dielectric^{50}. The enhancement factor values are obtained from FDTD simulations (We used commercial FDTD software of Lumerical Inc. for all FDTD simulations in this paper). A three-unit-cell-thick SAMM slab is used in the simulations. The enhancement factor is measured at 2 THz, which is a sufficiently low frequency and can be considered to be in the long-wavelength regime.

### Effective nonlinear electric susceptibilities calculation in Fig. 3

\(\chi _{{\mathrm{eff}},\,11}^{(3)}\) is retrieved from the simulations using a nonlinear transfer matrix method. The complex amplitude of the third-harmonic wave is obtained based on nonlinear FDTD simulations, whereas the linear electromagnetic property of the SAMM is obtained from separate linear FDTD simulations. A dispersion-less nonlinear dielectric (*ε*_{11} = *ε*_{33} = 6, \(\chi _{33}^{(3)}\) = 10^{−21} m^{2} V^{−2} for 10 GHz to 100 THz) and a lossy Drude model for gold^{49} are assumed as the constituent materials for Fig. 3a. Zinc oxide^{50} (*ε*_{11} = *ε*_{22} = 7.8, *ε*_{33} = 8.9, \(\chi _{33}^{(3)}\) = 10^{−18} m^{2} V^{−2} for 0.2 to 6 THz frequency) and the lossy Drude model for gold are assumed as the constituent materials, and a SAMM slab with two-unit cell thickness is used for the FDTD simulation for Fig. 3b. Zinc oxide^{51} (*ε*_{11} = *ε*_{22} = *ε*_{33} = 4, \(\chi _{33}^{(3)}\) = 10^{−18} m^{2} V^{−2} for 20 to 750 THz frequency) and aluminum^{52} are assumed as the constituent materials, and a SAMM slab with a half-unit cell thickness is used for the FDTD simulations for Fig. 3c.

### All-optical switch demonstration in Fig. 4

Owing to the limitation of the computational load, we obtained the results based on two-dimensional FDTD simulations (i.e. uniform in the *y*-direction in Fig. 4a). The pulse energy was calculated for a three-dimensional beam whose electric field profile in the radial direction is the same as that of the two-dimensional beam used in the simulations. The optical property of the GaN/AlN quantum well is calculated from two-level electron model^{11,34} whose carrier density, excitation energy, transition dipole length, relaxation time, dephasing time, and background relative permittivity are 5 × 10^{19} cm^{−3}, 0.8 eV, 0.5 nm, 100 fs, 10 fs, 5.37, respectively^{34,37}. Calculated intensity-dependent permittivity of the GaN/AlN quantum well is presented in Supplementary Fig. 4. The dynamics of the material property of the quantum well is treated as two-level electron model which incorporate equation of motion of electron population in sublevels^{53} in nonlinear FDTD simulations. Aluminum^{52} and silver^{52} are assumed to be the materials used for the metallic part of the SAMM and for the back-reflecting mirror, respectively. Silica is used for the upper and lower impedance-matching dielectric layers (relative permittivity: 2.074). A rigorous analytical model for the electric susceptibilities of SAMMs, a transfer matrix method, and a particle swam optimization are used for the structural parameter optimization.

### Modified Miller’s rule for the SAMMs in Fig. 5

The black dots in Fig. 5 are calculated using conventional Miller’s rule as

The red dots in the same figure are calculated using our modified Miller’s rule as

where

All \(\chi _{{\mathrm{eff}}}^{(1)}\) are obtained from the simulation results.

## Data availability

The authors declare that the data supporting the findings of this study are available within the article (and its Supplementary Information files) and are available from the corresponding author upon reasonable request.

## Code availability

The code to carry out simulations and analysis can be found at https://doi.org/10.6084/m9.figshare.12043308 and are available from the corresponding author upon reasonable request.

## References

- 1.
Kivshar, Y. S. & Agrawal, G.

*Optical Solitons: From Fibers to Photonic Crystals*(Academic Press, Cambridge, MA, 2003). - 2.
Chai, Z. et al. Ultrafast all-optical switching.

*Adv. Opt. Mater.***5**, 1600665 (2017). - 3.
Ou, J.-Y., Plum, E., Zhang, J. & Zheludev, N. I. Giant nonlinearity of an optically reconfigurable plasmonic metamaterial.

*Adv. Mater.***28**, 729–733 (2016). - 4.
Maruo, S., Nakamura, O. & Kawata, S. Three-dimensional microfabrication with two-photon-absorbed photopolymerization.

*Opt. Lett.***22**, 132–134 (1997). - 5.
Denk, W., Strickler, J. & Webb, W. Two-photon laser scanning fluorescence microscopy.

*Science***248**, 73–76 (1990). - 6.
Cai, W., Vasudev, A. P. & Brongersma, M. L. Electrically controlled nonlinear generation of light with plasmonics.

*Science***333**, 1720–1723 (2011). - 7.
Ren, M., Plum, E., Xu, J. & Zheludev, N. I. Giant nonlinear optical activity in a plasmonic metamaterial.

*Nat. Commun.***3**, 833 (2012). - 8.
Li, G. et al. Continuous control of the nonlinearity phase for harmonic generations.

*Nat. Mater.***14**, 607 (2015). - 9.
Tymchenko, M. et al. Gradient nonlinear pancharatnam-berry metasurfaces.

*Phys. Rev. Lett.***115**, 207403 (2015). - 10.
Chen, S. et al. Giant nonlinear optical activity of achiral origin in planar metasurfaces with quadratic and cubic nonlinearities.

*Adv. Mater.***28**, 2992–2999 (2016). - 11.
Boyd, R. W.

*Nonlinear Optics*(Academic press, Cambridge, MA, 2003). - 12.
Sirtori, C., Capasso, F., Sivco, D. L. & Cho, A. Y. Giant, triply resonant, third-order nonlinear susceptibility

*χ*^{(3)}_{3ω}in coupled quantum wells.*Phys. Rev. Lett.***68**, 1010–1013 (1992). - 13.
Lee, J. et al. Giant nonlinear response from plasmonic metasurfaces coupled to intersubband transitions.

*Nature***511**, 65–69 (2014). - 14.
Pendry, J. B., Holden, A. J., Robbins, D. & Stewart, W. Magnetism from conductors and enhanced nonlinear phenomena.

*IEEE Trans. Microw. Theory Tech.***47**, 2075–2084 (1999). - 15.
Zharov, A. A., Shadrivov, I. V. & Kivshar, Y. S. Nonlinear properties of left-handed metamaterials.

*Phys. Rev. Lett.***91**, 037401 (2003). - 16.
Kim, S. et al. High-harmonic generation by resonant plasmon field enhancement.

*Nature***453**, 757–760 (2008). - 17.
Rose, A., Larouche, S. & Smith, D. R. Quantitative study of the enhancement of bulk nonlinearities in metamaterials.

*Phys. Rev. A***84**, 053805 (2011). - 18.
Lassiter, J. B. et al. Third-harmonic generation enhancement by film-coupled plasmonic stripe resonators.

*ACS Photonics***1**, 1212–1217 (2014). - 19.
O’Brien, K. et al. Predicting nonlinear properties of metamaterials from the linear response.

*Nat. Mater.***14**, 379–383 (2015). - 20.
Jäger, I. & Fratzl, P. Mineralized collagen fibrils: a mechanical model with a staggered arrangement of mineral particles.

*Biophys. J.***79**, 1737–1746 (2000). - 21.
Gao, H. Application of fracture mechanics concepts to hierarchical biomechanics of bone and bone-like materials.

*Int. J. Fract.***138**, 101–137 (2006). - 22.
Smith, D. R. & Pendry, J. B. Homogenization of metamaterials by field averaging.

*J. Opt. Soc. Am. B***23**, 391–403 (2006). - 23.
Alù, A. First-principles homogenization theory for periodic metamaterials.

*Phys. Rev. B***84**, 075153 (2011). - 24.
Gorlach, M. A., Voytova, T. A., Lapine, M., Kivshar, Y. S. & Belov, P. A. Nonlocal homogenization for nonlinear metamaterials.

*Phys. Rev. B***93**, 165125 (2016). - 25.
Roke, S., Bonn, M. & Petukhov, A. V. Nonlinear optical scattering: the concept of effective susceptibility.

*Phys. Rev. B***70**, 115106 (2004). - 26.
Kleinman, D. A. Nonlinear dielectric polarization in optical media.

*Phys. Rev.***126**, 1977–1979 (1962). - 27.
Chang, T. et al. Broadband giant-refractive-index material based on mesoscopic space-filling curves.

*Nat. Commun.***7**, 12661 (2016). - 28.
Yang, X.-l & Xie, S.-w Expression of third-order effective nonlinear susceptibility for third-harmonic generation in crystals.

*Appl. Opt.***34**, 6130–6135 (1995). - 29.
Fredkin, D. R. & Mayergoyz, I. D. Resonant behavior of dielectric objects (electrostatic resonances).

*Phys. Rev. Lett.***91**, 253902 (2003). - 30.
Quinten, M.

*Optical Properties of Nanoparticle Systems: Mie and Beyond*(John Wiley & Sons, Hoboken, NJ, 2010). - 31.
Kim, J. Y. et al. Highly tunable refractive index visible-light metasurface from block copolymer self-assembly.

*Nat. Commun.***7**, 12911 (2016). - 32.
Chung, K., Kim, R., Chang, T. & Shin, J. Optical effective media with independent control of permittivity and permeability based on conductive particles.

*Appl. Phys. Lett.***109**, 021114 (2016). - 33.
Kim, R. et al. Metal nanoparticle array as a tunable refractive index material over broad visible and infrared wavelengths.

*ACS Photonics***5**, 1188–1195 (2018). - 34.
Suzuki, N. & Iizuka, N. Electron scattering rates in AlGaN/GaN quantum wells for 1.55-µm inter-subband transition.

*Jpn. J. Appl. Phys.***37**, L369 (1998). - 35.
Cho, S.-H. et al. A black metal-dielectric thin film for high-contrast displays.

*J. Korean Phys. Soc.***55**, 501–507 (2009). - 36.
Maniyara, R. A., Mkhitaryan, V. K., Chen, T. L., Ghosh, D. S. & Pruneri, V. An antireflection transparent conductor with ultralow optical loss (<2%) and electrical resistance (<6 Ω sq

^{−1}).*Nat. Commun.***7**, 13771 (2016). - 37.
Iizuka, N., Kaneko, K. & Suzuki, N. All-optical switch utilizing intersubband transition in GaN quantum wells.

*IEEE J. Quantum Electron.***42**, 765–771 (2006). - 38.
Hu, X., Jiang, P., Ding, C., Yang, H. & Gong, Q. Picosecond and low-power all-optical switching based on an organic photonic-bandgap microcavity.

*Nat. Photonics***2**, 185 (2008). - 39.
Koos, C. et al. All-optical high-speed signal processing with silicon–organic hybrid slot waveguides.

*Nat. Photonics***3**, 216 (2009). - 40.
Nozaki, K. et al. Sub-femtojoule all-optical switching using a photonic-crystal nanocavity.

*Nat. Photonics***4**, 477 (2010). - 41.
Martínez, A. et al. Ultrafast all-optical switching in a silicon-nanocrystal-based silicon slot waveguide at telecom wavelengths.

*Nano Lett.***10**, 1506–1511 (2010). - 42.
Ren, M. et al. Nanostructured plasmonic medium for terahertz bandwidth all-optical switching.

*Adv. Mater.***23**, 5540–5544 (2011). - 43.
Shcherbakov, M. R. et al. Ultrafast all-optical switching with magnetic resonances in nonlinear dielectric nanostructures.

*Nano Lett.***15**, 6985–6990 (2015). - 44.
Wood, R. M.

*Laser-Induced Damage of Optical Materials*(CRC Press, Cleveland, Ohio, 2003). - 45.
Ettoumi, W., Petit, Y., Kasparian, J. & Wolf, J. P. Generalized miller formulæ.

*Opt. Express***18**, 6613–6620 (2010). - 46.
Shin, J., Shen, J.-T. & Fan, S. Three-dimensional metamaterials with an ultrahigh effective refractive index over a broad bandwidth.

*Phys. Rev. Lett.***102**, 093903 (2009). - 47.
Choi, M. et al. A terahertz metamaterial with unnaturally high refractive index.

*Nature***470**, 369–373 (2011). - 48.
Munch, E. et al. Tough, bio-inspired hybrid materials.

*Science***322**, 1516–1520 (2008). - 49.
Ordal, M. A., Bell, R. J., Alexander, R. Jr, Long, L. & Querry, M. Optical properties of fourteen metals in the infrared and far infrared: Al, Co, Cu, Au, Fe, Pb, Mo, Ni, Pd, Pt, Ag, Ti, V, and W.

*Appl. Opt.***24**, 4493–4499 (1985). - 50.
Azad, A. K., Han, J. & Zhang, W. Terahertz dielectric properties of high-resistivity single-crystal ZnO.

*Appl. Phys. Lett.***88**, 021103 (2006). - 51.
Bond, W. L. Measurement of the refractive indices of several crystals.

*J. Appl. Phys.***36**, 1674–1677 (1965). - 52.
Palik, E. D.

*Handbook of Optical Constants of Solids*. Vol. 3 (Academic press, Cambridge, MA, 1998). - 53.
Taflove, A. & Hagness, S. C.

*Computational Electrodynamics: The Finite-Difference Time-Domain Method*. (Artech house, Norwood, MA, 2005).

## Acknowledgements

This work was supported by the Samsung Science & Technology Foundation (grant SSTF-BA1401-06) and an NRF (National Research Foundation of Korea) grant funded by the Korean Government (NRF-2017-Global Ph.D. Fellowship Program).

## Author information

### Affiliations

### Contributions

J.S. conceived the idea and supervised the project. T.C. and J.S. conducted the theoretical analyses. T.C. performed the numerical simulations. T.C., S.J., and M.H. performed the matched dielectric layer design.

### Corresponding author

## Ethics declarations

### Competing interests

The authors declare no competing interests.

## Additional information

**Publisher’s note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

## Supplementary information

## Rights and permissions

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

## About this article

### Cite this article

Chang, T., Jeon, S., Heo, M. *et al.* Mimicking bio-mechanical principles in photonic metamaterials for giant broadband nonlinearity.
*Commun Phys* **3, **79 (2020). https://doi.org/10.1038/s42005-020-0352-0

Received:

Accepted:

Published:

## Comments

By submitting a comment you agree to abide by our Terms and Community Guidelines. If you find something abusive or that does not comply with our terms or guidelines please flag it as inappropriate.