Simplified model for periodic nanoantennae: linear model and inverse design

We determine and use a minimal set of numerical simulations to create a simplified model for the spectral response of nanoantennae with respect to their geometric and modeling parameters. The simplified model is then used to rapidly obtain best-fit modeling parameters to match experimental results, accurately predict the spectral response for various geometries, and inversely design antennae to have a desired performance. This method is structure and model independent, and is applied here to both nanoantenna pair arrays and strips modeled using a 3D finite-element method and 2D spatial harmonic analysis, respectively. Typical numerical simulations may need hours per model, whereas this method, after the initial time to obtain a baseline set of simulations, requires only seconds to analyze and generate spectra for new geometries.

The performances of fabricated nanoantennae typically differ from ideal models due to effects not accounted for in the model. The actual size and shape of the antennas differ from the initial design due to fabrication limitations and systemic error. Additionally, the properties of the metal may differ from those used in simulations. For example, the permittivity of gold is typically taken from bulk film measurements, but does not account for gold lattice defects, chemical size effects, and most importantly, surface roughness, which can have a significant effect on nanoparticle resonances [28]. Therefore, in order to match an experimental result with a simulation, the parameters of the simulation must be changed to 'effective' values which account for these discrepancies.
[24] Many time consuming simulations may be required to find parameters, which match the experimental result.
Additionally, pre-fabrication simulations are required to determine which dimensions will lead to a desired performance. Inversely designing a nanoantenna using nature-inspired optimization methods [29] can do this, although feedback from fabrication is necessary in order to account for non-ideal fabrication constraints and variations in material properties as explained above, while also being time consuming, especially for a 3D simulation.
In order to reduce the time and computer resources required, we create a simplified equivalent model of a nanoantennae system, which results in the ability to analyze many geometries with minimized computational effort. We apply an equivalent model to both nanoantenna pair arrays and strips using an experimental design methodology [30][31][32][33][34]. This results in a model based on simple linear equations, which accurately describe the antenna system over a limited parameter range, and although illustrated here with paired nanoparticles and strips, may be generally applied to any simulation.

Nanoantennae fabrication and modeling
The paired nanoantennae diagramed in Fig. 1 (first published as Sample A in [24], but herein identified as Sample 1) were used in addition to the nanostrips shown in Fig. 2. The nanostrips were fabricated on a glass substrate coated with 15 nm of indium-tin-oxide (ITO). E-beam lithography was used with a PMMA resist to pattern the nanoantennae. The substrate was coated with 5 nm of titanium as an adhesion layer followed by the periodic gold nanoantennae array. In Sample 1, the nanoantenna pair array is the array of paired elliptic-shape gold particles, and sample 2, nanoantenna strips, consists of a periodic array of paired gold strips. Both of these structures are resonant in optical frequencies. Transmission and reflection spectra over the visible range are taken using linearly polarized light; incident polarization across the gap induces the primary resonance (x direction; P polarization), and polarization parallel to the gap induces a secondary resonance in the case of pair arrays, or no resonance in strips (ŷ direction; N polarization).  Typically, the post-fabrication geometry is found from SEM measurements with limited resolution. More precise effective dimensions are determined by matching modeled transmission and reflection spectra with the experimental results. This 'fabrication-calibrated' model, using the proper effective dimensions, can then be used to retrieve optical properties such as index of refraction, and the magnitude and wavelength dependence of the local field enhancement.

Simplified system
We create a set of simulations, which is used to determine a simplified system describing the results. In our case, a linear equation with variable interactions is used, which once defined, may then quickly predict the results of future simulations. The simulation has n input parameters 1 ,..., n P p p   =   (these could be model variables such as width, height, etc.), output f , and coefficients β relating the two. A general polynomial system is shown (up to quadratic terms) in Eq. (1); more complex higher-order terms may be added if necessary.
The methodology we follow is a well developed approach called experimental design [34]. This method is widely used to analyze systems with multivariate input using a limited number of experiments, and is typically applied when experiments are costly, such as in biological research, and although the method has existed for some time, its application to metamaterial electromagnetic modeling has not been explored. Herein we use a two-level fractionalfactorial [30-33] experimental design.
In order to use this system to predict the output f , we must first determine β . The solution to this system may be more readily understood using matrix notation. For the linear plus interaction system considered here, there is one constant term, n linear terms, and interaction terms for a total of ( ) (a set of t input parameters for each of m simulations). β is a t -dimensional vector of coefficients (a coefficient for each term), and the output f is an m -dimensional vector (a scalar result for each of the m simulations). Note that although mathematically the system input is X , it is determined from the fundamental model parameters P , e. g. the m th row of X is given by: We use a two-level design for P , so that over the m simulations each input parameter i p is either a low or high value, -1 or 1 (e.g. width is 90 nm or 110 nm). In order to get all t terms of vector β , we must run enough simulations so that each of the t columns of X are linearly independent, i.e., the system is not under-determined. For instance, if we have 5 input parameters, then 16 t = and we must design 16 sets simulation parameters m P , 1 16 m = … , so that we have 16 linearly independent simulations with X m , 1 16 m = … , and thus, X , is a 16 16 × matrix. This is easily done with common statistical programs, or for example with Matlab's 'fracfact' function [35].With m P properly defined (which in turn gives m Χ ) , we then run the simulations to retrieve the results m f , so we may finally solve for all coefficients β using β = X f . With vector β obtained, the simplified system is determined. From now on, we may use Eq. (2), for any set of input parameters P (which therefore defines Χ , a 1 t × matrix), to predict the output f (a scalar value).

Simplified nanoantennae model
Applying experimental design methodology to nanoantennae modeling reduces the number of full simulations that must be run to a reduced set, which may then be used to describe the system. Any number of input variables n may be used; potential input variables for the current samples include geometric parameters such as gap, width, height, and periodicity, or other parameters such as loss factor. Our initial goal is to find the model parameters which result in spectra that best match the experimental data. In this work, we use 3D FEM and 2D SHA simulations, however it is important to note that since the analysis only relies on input parameters P versus an output value f , the simulation itself is a 'black box' and may consist of any modeling technique desired. While f may be defined as one of many potential output values, we define it using one of two separate but complimentary methods.

Fit parameters for nanoantennae
In the first method, each set of simulated spectra (transmission and reflection at both polarizations) are defined according to five output values. These values are specifically chosen to describe nanoantennae spectra, and may not apply to other types of models. (error in full-width-half-max for P polarization). Since the resonance position error is on the order of tens of nanometers, but magnitude error is between 0 and 1, the five output values are normalized to a range from -1 to 1. We therefore have five f output values for each simulation. A separate simplified model (separate set of coefficients β ) is determined independently for each of the five fit criteria, meaning that we can predict each fit parameter for a given input P . In order to achieve the overall best fit between the experiment and a simulation, the input values P are found to minimize a weighted root-mean-square of the five fit criteria, F , Eq. (3). A weight is ascribed to each of the fit parameters to allow for a difference in the importance of each fitting criterion. Minimizing Eq.
Using this simplified system of five fit parameters also, and perhaps more importantly, allows us to reverse the process and inversely design the nanoantennae. Instead of providing experimental data to determine each simulation's error, 1 f to 5 f , we may instead provide desired values for an optimal nanoantenna, (e.g. a desired resonance wavelength). In this case minimizing Eq. (3) determines the values for the model, P , which result in the desired performance. Using this method, many antennae may be designed rapidly and the results explored, simply through minimizing a simple equation. This results in effective and accurate generation of design parameters for a desired antenna performance.

Fit per wavelength
A second method is to define an individual output f as the percent transmission or reflection for each spectrum at each wavelength. For example, suppose each simulation outputs four spectra (transmission and reflection for P and N), and wavelengths of 400 nm to 900 nm in steps of 20 nm (25 steps), then we would have 100 different f values, and therefore 100 β vectors. This complete set of coefficients allow us to recreate or predict the full spectra, comparable to an FEM or SHA model, for any set of input parameters P . However, since these spectra are the result of solving simple equations, the spectra are obtained nearinstantaneously. This method, similar to using the 'fit parameters' method, may also be used to fit experimental data. Model input values P may be determined which minimizes the rootmean square error between the entire experimental spectra and the predicted spectra. This method has the advantage over the 'fit parameters' method that it is not model specific, and will fit spectra for any model, however it cannot be used for inverse design.

Inverse design
A combination of both methods leads to a complete inverse design system. The performance of the antenna (primary and secondary resonance wavelengths, resonance magnitudes, etc.) may be manually defined, and the 'fit parameters' method may be used to determine the optimal model values, P , to achieve this result. Similar to more complex inverse design optimization methods [29], our simple model may be used to find an optimal design, and additionally if the simplified model were built from an array of experimental samples, instead of simulations, the fabrication and material variations that limit stochastic optimization methods would be built inherently into our simplified system, whether or not the deviations were measured, or even understood. After determining an optimized geometry, the 'per wavelength' method may then be employed to generate an entire set of predicted spectra for the inversely designed sample.

Nanoantennae fitting
The nanoantenna pair array (sample 1)[24] as described in Fig. 1, has    The ability of the 'fit per wavelength' method to predict entire spectra was used to generate the spectra for the given best-fit values, and is shown in Fig. 3 compared to the experimental data. In order to validate the ability to predict an entire set of spectra using method 2, the predicted best-fit curves are also compared to actual the FEM models, which were simulated using the above best-fit values. The maximum error between the actual FEM model and the simple linear model prediction spectra averages less than 1%, with a peak of 4% near the resonance. This indicates that a linear model is accurate over the range of parameters used, and agrees well with the predicted spectra.
In addition to achieving a best-fit model to experimental data, the effects of each variable on the response of the antenna may be analyzed. For instance, if we consider the slope of each term of 1 f ( β coefficients in Table 1), the function describing the error in the wavelength of the primary resonance, with respect to each input variable, we may understand which dimensions are essential to the robustness of the design, and to the response of the nanoantennae. We immediately see that, in order of importance, the long-axis and the gap have the greatest effect on the primary (P) resonance location, whereas 2 f , the resonance magnitude is dominated by the loss factor, α , as might be expected. Fig. 3. Nanoantenna array (sample 1), transmission and reflection for primary polarization (P) and secondary polarization (N). Experimental results vs. method 1 (fit parameters) best fit. Predicted spectra were generated using method 2 (fit per wavelength). The results of the actual 3D FEM model are shown using the same best-fit parameters

Nanoantennae inverse design
Instead of using the existing fabricated sample, we also inversely design an optimum nanoantennae structure by manually defining the values , , , ,  [12,118,55] nm. The predicted spectra show a strong resonance with only 10% transmission at 690 nm as shown in Fig. 4, and the result is closely matched by an FEM model with the same parameters.

Nanoantenna strips fitting
A similar analysis was performed on sample 2, nanoantenna strips, as shown in Fig. 2. Fig. 5. Nanoantenna Strips (sample 2), transmission and reflection for primary polarization (P) and secondary polarization (N). Experimental results vs. method 1 (fit parameters) best fit. Predicted spectra were generated using method 2 (linear fit per wavelength). The results of the actual 2D SHA model are shown using the same best-fit parameters, and agree well with the predicted curve.
From the QCM and SEM image, the dimensions [ , , , , ]