Parametric study of optical transmission through plasmonic hole arrays modulated by the phase transition of vanadium dioxide

This Supplementary Document includes additional experimental and computational details on the fabrication, measurements and simulations of the noble-metal–vanadium-dioxide bilayer arrays of periodic nanoholes, as well as five of the figures referenced in the main text.


Introduction
The discovery of the extraordinary optical transmission (EOT) phenomenon in 1998 [1] gave impetus to the field of nanoplasmonics [2][3][4][5], which deals with confining and manipulating light in the vicinity of metallic nanostructures [6][7][8] and has generated a myriad of applications [9][10][11][12][13][14]. To realize the full potential of nanoplasmonic structures, it is necessary to dynamically control their optical responses [15][16][17][18][19][20][21]. By varying periodicity, film thickness and hole diameter, we demonstrate an iterative strategy for maximizing the EOT switching induced by the semiconductor-to-metal phase transition (SMPT) of vanadium dioxide (VO 2 ) in Au+VO 2 and Ag+VO 2 bilayer hole arrays. We also trace the origin of the 'reverse' sense of the switching mechanism (i.e., lower T 00 in the semiconducting phase) to light-trapping Fabry-Perot resonances being excited more efficiently in the semiconducting phase of the VO 2 film due to its nearly-flat optical constants within the explored spectral range.

Experimental measurements
For the bilayers, sub-stoichiometric VO x films were deposited by pulsed laser deposition (PLD) on quartz substrates by 266-nm laser pulses with a repetition rate of 10 Hz, incident onto a 99.99%-purity vanadium (V) metal target (from Alfa Aesar) under a 5-mTorr partial O 2 pressure at room temperature. The samples were then annealed in a furnace for 120 min under 280 mTorr of O 2 at 450 °C to crystalize the amorphous VO x into polycrystalline VO 2 films. For the VO 2 film in Fig. 1(a), the PLD wavelength, repetition rate, V-target and annealing O 2 pressure were 248 nm, 15 Hz, 99.9%-purity (from PLDtargets.com) and 250 mTorr. Following the VO 2 synthesis process, optical transmission was measured as a function of varying temperature in order to confirm the SMPT of the VO 2 films. The desired VO 2 thickness was confirmed by measuring the step height by optical (Filmetrics Profilm3D) or stylus (KLA-Tencor P-6) profilometry. A 200-nm-thick Au layer was deposited on top of the VO 2 layer via RF sputtering (Perkin-Elmer 2400, 13.6 MHz) at a base pressure of 6.5×10 -6 Torr, 100-watt RF power, 15-sccm Ar flow rate and deposition rate of 30 nm/min. Subsequently, 100×100 μm 2 arrays of holes of different periodicities were milled through the bilayer using a single-column focused ion beam (Hitachi FB-2000A FIB). The periodicity and hole diameter were specified in the Nanometer Pattern Generation System (JC Nabity Lithography Systems), and the beam-current dose was adjusted by the dwell time.
The optical setup for measuring the zero-order transmission (T 00 ) is shown schematically in Fig. S1(f, inset). Once the incident light was focused on the nanohole pattern using an xyzstage, the beam spot size was adjusted via a pinhole to ensure all of the illumination fit within the pattern. Another lens was placed on the output side of the sample to collect and send the transmitted light through a multimode optical fiber to a grating spectrometer (Andor 193i) equipped with Si and InGaAs detectors. Performing grain analysis in software (Gwyddion v2.53) on one of the higher-magnification SEM images, we estimated the average openingaperture diameter of holes in the Au layer to be D Au(-holes) ≈ 230 nm [ Fig. S1(b)], while the corresponding apertures in the VO 2 layer came out significantly narrower with an average D VO2(-holes) ≈ 125 nm [ Fig. S1(c)]. The average periodicity in the xy-plane (i.e., sample surface) was determined by computing (in Gwyddion) a two-dimensional fast-Fourier transform (2D-FFT) [ Fig. S1(d)] of the 10k× micrograph in Fig. S1(a) and analyzing the extracted reciprocal-space line profiles in the x-and y-directions [ Fig. S1(e)]. Both measured and simulated T 00 spectra for the 650-nm-period bilayer hole array exhibit reverse switching, i.e., higher peak transmission for MetVO 2 than SemiVO 2 , in contrast to transmission through an unperforated VO 2 film [cf. Fig. 1(b)].

Numerical simulations and Fano model
The computational results presented in this paper derive from a subset of some 1000 simulations spanning about 6000 runtime hours. These totals include many variations not shown here, such as most of the results for Ag+VO 2 , Au+VO 2 without a substrate, preliminary runs and convergence tests. Most of the simulations were performed on a single workstation (Dell Precision T5600) with two eight-core processors (Intel Xeon E5-2650, 2.0/2.8 GHz base/turbo frequency) and 64 GB RAM, running Windows 10 Professional (64bit).
We chose the finite-difference time-domain (FDTD) method for some of the simulations [e.g., Fig. S1(g, h) and Fig. 2(f)] because the time-domain algorithm executes a broadband wavelength sweep in one run, which takes approximately the same amount of time (~12 hours) as a single-wavelength simulation, whereas the frequency-domain finite-element method (FEM) would require a separate run for each desired wavelength. On the other hand, this decoupling of spectral points gives FEM an enormous time-saving advantage when only a few wavelengths are needed for each parameter value to find the maximum transmission (peak-T 00 ). If, for example, all 304 simulations (152 film thickness values in each VO 2 phase) that went into generating Fig. 5 had been performed with FDTD instead of FEM, the runtime hours would have been more than 3500 instead of 20.

Finite-difference time-domain (FDTD) method
The dispersive constitutive relations, which in the frequency domain are products of harmonic functions representing materials and fields, become convolution integrals in the time domain; therefore, analytical models of frequency-dependent permittivities (and, if applicable, permeabilities) must be implemented. Once the fields have been explicitly solved for in an alternating ("leapfrog") fashion by a time-stepping algorithm, the frequency-domain responses can be obtained via Fourier transforms.
The simulation domain represents a periodic bilayer hole array of infinite extent in the xyplane, sandwiched between a semi-infinite air (n air = 1.00) medium on the input side and a semi-infinite glass (n glass = 1.50) or air medium on the output side. The incident illumination is an x-polarized plane wave launched as a broadband pulse of 7.98-fs duration and 549.62-THz bandwidth (400-1500 nm). A field monitor is placed in an xy-plane 50 nm away from the Au-air interface into the air medium and another one 50 nm away from the VO 2 -glass interface into the substrate; these monitors are used for calculating, by means of Lumerical's built-in grating and transmission script functions, the fractions of the incident power reflected and transmitted (T 00 ) into the zeroth diffraction orders. Simulation time is set to at least 3000 fs to allow resonant responses to dissipate after the incident pulse interacts with the system.

Finite-element method (FEM)
The FEM we employ solves Maxwell's equations by the variational approach of minimizing an integral functional of the field variables with respect to small variations in these variables. In equilibrium, the physical quantity of interest (e.g., resonant frequency or transmitted power) can be expressed as a parameter that produces a stationary value of the minimized functional [22]. Inserting the continuity-matched local interpolation functions into the variational equation yields a global matrix system of equations whose solution is a piecewise polynomial approximation of the exact solution [23].
Working within the Wave Optics module of COMSOL Multiphysics, we again take advantage of the geometrical symmetries of the xy-periodic hole array and infinite-extent illumination source (x-polarized plane wave at normal incidence in the z-direction), which allow us to cut the runtime by simulating only a quarter of the structure [e.g., see Fig. 12(c)] without loss of information. The perfect electric conductor walls at the yz-boundaries and perfect magnetic conductor walls at the xz-boundaries render the ¼-structure equivalent to an infinite xy-periodic array of holes. The incident light is injected from a wave-on port boundary placed ~500 nm away from the air-Au interface within the air subdomain. At a wave-off output-port boundary within the substrate subdomain, ~500 nm away from the VO 2glass interface, the zero-order power transmission is calculated from the S-parameters of the scattering matrix as T 00 = |S 21 | 2 , by means of the COMSOL expression abs(ewfd.S21)^2.

Thickness of VO 2 layer
Sections 4.2 and 4.3 revealed the VO 2 -layer thickness to be one of the key geometrical parameters affecting the EOT modulation. Therefore, we compare below the intensity images and power flow vectors for two specific t VO2 values, 494 nm [ Fig. S4(a, b)] and 700 nm [ Fig.  S4(c, d)], while keeping the other parameters fixed. Fig. S4. FEM images of electric-field intensity in dB, 20 log 10 (E/E 0 ), where Ε 0 is incident electric field at input port, with superimposed arrows of power flow S in W/m 2 , for Au+VO 2 hole arrays on glass with different VO 2 -layer thickness (t VO2 ). Array period (720 nm), thru-hole diameter (290 nm) and Au-layer thickness (200 nm) are identical for both cases. Illumination source is x-polarized plane wave (E inc || î) injected from air medium (not shown) at normal incidence in upward z-direction. Color scale ranges from dark blue (+30 dB) to yellow (+90 dB) to dark red (+130 dB) and represents intensity in xz-and yz-planes of symmetry bounding simulated ¼-cylinder hole (see These values are selected because, as seen in Fig. 5(a, bottom panel), the two corresponding Air-Au+VO 2 -Glass hole arrays exhibit comparable MetVO 2 peak-T 00 values (8×10 -5 vs. 4×10 -5 ) but 100-fold-different SemiVO 2 peak-T 00 values (0.08×10 -5 vs. 10×10 -5 ). In particular, the SemiVO 2 transmission dips dramatically and M2S-ratio = 101 at λ peak = 773 nm [see Fig. 5(a, d)] for t VO2 = 494 nm, which we attributed earlier to the excitation of a vertical FP-type anti-resonance. When the thickness is increased to t VO2 = 700 nm, M2S-ratio = 0.37 at λ peak = 773 nm [see Fig. 5(a, e)]. Despite the large difference in M2S-ratios for these two t VO2 values, the two sets of intensity images in Fig. S4 look surprisingly alike pairwise, both in the xz-plane [(a) vs. (c)] and yz-plane [(b) vs. (d)], and especially in the Au, VO 2 and hole media. Yet, a possible sign of the FP anti-resonance responsible for the SemiVO 2 dip at t VO2 = 494 nm in Fig. 5(a, bottom panel) may be the narrow dark-blue feature located just above the VO 2 -Glass interface in Fig. S4(b, left panel). [A seemingly equivalent low-intensity spot in Fig. S4(d, left panel) actually occurs farther to the left of the aperture rim and higher above the VO 2 -Glass interface.]

Array period with scaled thru-hole diameter
Finally, let us briefly return to the first optimization sweep (Section 4.1), varying P array (= 3D thru-hole ), in order to compare the electric-field intensity and power flow for three representative Air-Au+VO 2 -Glass hole arrays whose M2S-ratios straddle the Fano-type curves in Fig For P array = 720 nm and D thru-hole = 240 nm, M2S-ratio = 8.8 at λ peak = 759 nm [see Fig.  4(b) and Fig. 3(a, third panel from left)] and, as seen in Fig. S5(c, d), the overall intensity and power flow in the glass medium above the MetVO 2 layer (right panels) clearly exceed those above the SemiVO 2 layer (left panels). The power flow also appears to spread out laterally more within the SemiVO 2 layer than within the MetVO 2 layer, which presumably results in less light being available for T 00 .
For P array = 1110 nm and D thru-hole = 370 nm, M2S-ratio = 2.8 at λ peak = 1137 nm [see Fig.  4(b) and Fig. 3(a, rightmost panel)] and, as seen in Fig. S5(e, f), the overall intensity and power flow in the glass medium are higher than in the other two cases, which, as mentioned in Section 3.3, comes about because the EOT peak is redshifted towards the near-IR where Au is a better plasmonic metal. However, it is somewhat counterintuitive that, at least visually, the light penetrates much deeper into the SemiVO 2 layer [left panels in Fig. S5(e, f)] than it does in the MetVO 2 layer [right panels in Fig. S5(e, f)], and yet M2S-ratio ≈ 3 only. In the preceding case [ Fig. S5(c, d)], the visual differences in how much light penetrates the VO 2 layer are less obvious but M2S-ratio ≈ 10. It seems that, when the holes are large enough to funnel the bulk of the transmitted light, the lateral penetration of the fields into the VO 2 film becomes less of a factor in the Met-to-Semi switching ratio. The difficulty of directly correlating the visual comparison of electric-field intensity (or Poynting vector) in each VO 2 phase with the numerical comparison of peak-T 00 in each phase serves as another reminder that the mechanism of the SMPT-induced EOT modulation is in fact richer and subtler than it appears. Fig. S5. FEM images of electric-field intensity in dB, 20 log 10 (E/E 0 ), where Ε 0 is incident electric field at input port, with superimposed arrows of power flow S in W/m 2 , for Au+VO 2 hole arrays on glass with three different periods (P array ) and scaled thru-hole diameters (D thru-hole = P array /3). Au-layer thickness (200 nm) and VO 2 -layer thickness (200 nm) are identical for the three cases. Illumination source is x-polarized plane wave (E inc || î) injected from air medium (not shown) at normal incidence in upward z-direction. Color scale ranges from dark blue (+30 dB) to yellow (+90 dB) to dark red (+130 dB) and represents intensity in xz-and yz-planes of symmetry bounding simulated ¼-cylinder hole (see Fig. 12). Left panel of each subfigure shows half of unit cell cross section in SemiVO 2 phase; right panel mirrors same geometry but in MetVO 2 phase. Plots for each pair of P array and D thru-hole values in (xz, yz)-planes [see also Fig. 3 and Fig. 4]: (a, b) 510 nm and 170 nm, at λ peak = 590 nm; (c, d) 720 nm and 240 nm, at λ peak = 759 nm; and (e, f) 1110 nm and 370 nm, at λ peak = 1137 nm.