In the last few years, the complexity of systems focused on materials science and experimental physics is gradually increasing to the point that analytical approaches to theoretical modeling are being overshadowed in the field by their numerical counterparts. This trend in research is thus promoting the evolution of comprehensive numerical packages that provide extremely accurate descriptions of physical processes.[1-3] Handily, commercial packages now provide highly reliable modeling of phenomena like structural or fluid behaviors,[4-7] wave propagation,[8-13] and thermal transport,[14-16] to name but a few. As a consequence, the availability of such powerful tools is modifying the way research is performed, eventually suggesting new possibilities to optimize the scientific procedure itself. Within the portfolio of experimental characterization techniques, it is remarkable what detailed insight can be acquired when performing spectroscopic ellipsometry of thin films and bulk materials,[17-19] including the monitoring of film growth,[20] surface roughness,[21] detection of micro- and nano-structures on a surface,[22-24] and determination of general optical material characteristics.[25-28] Indeed, after design and realization of a specific system, ellipsometric characterization typically represents a crucial test-bed, providing a validation of the expected functionalities or helpful hints for improvement. At this stage, a trial-and-error procedure to achieve the desired target very often results in numerous attempts. Considering the vast range of applicability of ellipsometry, the possibility of implementing this technique as a reliable numerical method, able to a priori predict the result of the real experimental analysis before actually performing it, is evidently extremely convenient. This is especially true when the focus of the study is a complex nanostructured system whose fabrication can be largely time-consuming and expensive.In literature, preliminary attempts at implementing a numerical ellipsometric analysis by means of finite element method (FEM) simulations in different spectral ranges have been reported,[29-33] mainly proposing an ideal starting point such as solution of single layer materials with functionalities in the VIS-NIR range. However, the approach in these cases is limited to a scattered field study instead of a full field analysis.Here, we present a novel method, implemented using the COMSOL Multiphysics FEM commercial platform,[34] providing a comprehensive ellipsometric analysis of general multi-layer systems with extreme freedom of design in terms of thickness, composition and number of layers. For each considered system, the proposed numerical ellipsometric analysis (NEA) allows for the calculation of the main physical quantities of interest for different polarizations and angles of incidence.
The numerical simulation environment
In the simulation geometry presented in this work (Fig. 1a), an electromagnetic plane wave (E(r̄,t) = E0 exp[i(k̄·r̄ − ωt + ϕ)]) impinges on a generic multi-layer system from the superstrate (typically air), with a specific incident angle and polarization (s-polarization (TE) or p-polarization (TM)). Light interacting with the structure is then collected through both its substrate (typically glass) and superstrate, permitting us to compute the optical quantities of interest (transmittance (T), reflectance (R), the Brewster angle, field maps and so on). This is done by solving the frequency-domain partial differential equation (PDE) that governs the E and H fields associated with the electromagnetic wave propagating through the structure:
Fig. 1
(a) Pictorial view of a basic ellipsometer setup. The sample is a thin film on a glass substrate. The sketches of the electric and magnetic field components for s- and p-polarized light impinging on the sample are depicted in the case of a 2D (b) and a 3D (c) numerical simulation.
With the conditions of the electric conductivity σ = 0 and non-magnetic materials (μr = 1), the previous equation reduces to:Here k0 is the incident wavevector in a vacuum or in air (k0 = 2π/λ), while εr represents the dielectric permittivity of the material. By providing as an input, the values of real and imaginary parts of the refractive index of any considered material, the software retrieves the corresponding dielectric permittivity (εr) and numerically solves eqn (2) to obtain the E field distribution. Later on, from the field values, the relevant quantities can be derived to numerically implement the NEA. Indeed, in the case of a standard ellipsometric analysis, we have to measure the ratio ρ = r̃p/r̃s, where r̃p and r̃s are the complex Fresnel coefficients related to the square roots of the reflectances for the two different polarizations, and . Asρ is a complex quantity, it can also be written as tan(Ψ)ei. Thus, tan(Ψ) is related to the amplitude ratio upon reflection (Ψ = arctan(rp/rs)), whereas Δ represents the phase shift between the two components (Δ = ϕp − ϕs).[33,35] Since ellipsometry measures the ratio (or difference) of two values (rather than their absolute values), it is very robust, accurate and reproducible.In order to simulate the behaviors of TE (s-pol) and TM (p-pol) polarized light[36] in a 2D environment, it is necessary to properly write the components of the electromagnetic fields with respect to the incidence plane (xy). To select the s-polarization in COMSOL, an “out of plane” configuration has to be used (E = (0, 0, 1)), while to select the p-polarization, the user has to apply the “in plane wave-vector” scheme; in this case the magnetic field has to be used as H = (0, 0, 1). In a 3D simulation, the incident plane corresponds to the xz plane and the polarizations are set accordingly: E = (0, 1, 0) for s-pol or H = (0, 1, 0) for p-pol. The four sketches shown in Fig. 1b and c clarify how to set the field components in COMSOL.To optimize the computation time of the NEA calculation, Floquet periodic boundary conditions (PBCs) were set to all boundaries. The PBCs ensure the simulation of systems with infinite size along the selected directions (ESI, Fig. S1†). Another way to optimize the computation process consists of parallelizing the solvers by introducing two “frequency domain studies”, able to solve eqn (2) for both s- and p-polarized light.The discretization of the problem was realized through the choice of a dense mesh weave, ensuring that the numerical simulation takes into account the size of the real layers (sub-nanometric sizes are not considered to avoid unphysical results, Fig. S1c†).
Results and discussion
In order to validate the effectiveness of the NEA modeling, numerical predictions were compared with the results of measurements performed on corresponding real samples. Several multi-layer systems have been designed, fabricated and characterized using experimental ellipsometric analysis. The first and simplest system was a single Ag layer deposited on a glass substrate (20 nm thickness). The following system was a hyperbolic metamaterial[37] made of a stack of five alternated Ag/ITO bi-layers (ITO stands for indium tin oxide) with the single bilayer having a thickness of 40 nm (20/20). The third system was characterized by a three-materials unit cell repeated three times. The unit cell was made of Al2O3 doped zinc oxide (AZO), ITO and Ag. The numerical task was particularly difficult because each of the three cells that were the constituents of the system had different thicknesses for each layer. The last and most complex system comprised an asymmetric optical cavity, with nanoscale features, also known as metal–insulator–metal (MIM).[38,39] Starting from the glass substrate, this system was composed of an Ag layer covered by a thick ITO slab, representing the dielectric cavity, and an Ag/thin-ITO bilayer on top. All samples were experimentally fabricated by exploiting a DC sputtering technique. Then, for all cases, the experimental reflectance, transmittance and ellipsometric angles Ψ and Δ were measured by means of a M-2000 ellipsometer (J.A. Woollam).Fig. 2a and b respectively show the directly measured reflectance and transmittance, as well as the corresponding numerical curves provided by the NEA simulation for the first case study. The reflectances, Rp and Rs curves, were measured and calculated by considering an incident angle of θi = 50°, while the transmittance was measured and simulated at normal incidence. The minimum value at about 350 nm in the reflectance curves (and the related maximum in the transmittance ones) is referred to as the Ferrell–Berreman mode for a silver nanometric layer.[40] The difference of about 20% between the amplitude of the experimental and simulated curves is due to the presence of the glass substrate in the real sample, that is partially considered in the simulated case. Indeed, in order to optimize the computational time, the size of the glass in the simulation was decreased to 1200 nm instead of the real 1.1 mm. The results are different only in absolute values but not in the general trends of the curves. Apart from that, the agreement between measurements and simulations is quite satisfactory. Fig. 2c and d show the Ψ and Δ behaviors, respectively.
Fig. 2
Reflectances (Rp, and Rs, red and black lines with symbols) and transmittance (T, blue lines with squares) measured (a) and calculated using the NEA model (b) for a single Ag layer on a glass substrate. Comparison between Ψ and Δ for the experimental (colored triangles and squares) and numerical (black lines) cases respectively, for the same samples (c and d).
The second test of our NEA model was conducted on a hyperbolic metamaterial (HMM) that, by alternating lossy metal layers with dielectric ones, acquires particular optical features. Thickness and composition of the nanolayers can be opportunely designed in order to exploit unusual optical properties in a desired spectral range. In fact, the particular design and choice of sizes and constitutive materials leads to the epsilon-near-zero-and-pole (εNZP) HMM,[41-43] allowing extraordinary light confinement properties. The proposed HMM system is sketched in the inset of Fig. 3b and is composed of five Ag/ITO bilayers, characterized by a fill fraction of 50%. Fig. 3a and b show, respectively, the experimental and numerical reflectances (Rp, and Rs) and the transmittance (T) calculated by the NEA model together with the comparison between experimental and numerical Ψ and Δ for the same sample (Fig. 3c and d). Also in this case, all measured quantities show impressive agreement with their numerical counterparts. Fig. 3e and f show the electric field maps respectively for an s-pol and a p-pol wave at λ = 390 nm, where it is clearly evident that the p-pol wave is able to penetrate more efficiently through the HMM structure than in the s-pol case, where the wave is transmitted by the medium in accordance with the metamaterial prediction.[44-47]
Fig. 3
Reflectances (Rp, and Rs, red and black lines with symbols) and transmittance (T, blue lines with squares) measured (a) and calculated using the NEA model (b) for the five bilayers Ag/ITO HMM, as sketched in the inset. Comparison between experimental (colored lines with symbols) and numerical systems (solid lines) Ψ (c) and Δ (d), for the same sample. (e and f) Electric field maps distribution for TE (s-pol) and TM (p-pol) waves respectively, extracted for an impinging wavelength of 390 nm.
A further and more complex test of the NEA model has been performed by considering the third system where every kind of thickness periodicity is absent. This choice reflects a more typical experimental situation with different thicknesses of all the layers and materials involved. The sample was composed of a unit cell of three materials, Ag, ITO and AZO, repeated three times, where the thickness of each layer was different (see inset of Fig. 4b). Starting from the glass substrate, the layer thicknesses were: 19 nm (Ag), 15 nm (ITO), 15 nm (AZO), 20 nm (Ag), 21 nm (ITO), 16 nm (AZO), 23 nm (Ag), 27 nm (ITO) and 15 nm (AZO). Fig. 4 reports the results related to this sample; the agreement between experiments and simulations is also quite remarkable in this case. The numerical ellipsometer simulation represents a real chance to overcome the typical theoretical approach based on the Transfer Matrix Method (TMM), that becomes particularly cumbersome when an additional 2D or 3D periodical structure is placed on top of the most external metamaterial surface. With this aim, the user-friendly interface allows positioning, without particular effort, nanostructures of different shapes and sizes on the top surface, and the proposed NEA tool allows the study of how the transmitted signal is affected by them.[48-51] In the fourth considered case, we simulated the optical behavior of an asymmetric multi-layer that does not respect the typical metamaterial fill fraction condition. In fact, when one of the layers is much thicker than the other ones, nanometric cavity behavior comes in, enabling confinement and re-direction of light. This kind of metamaterial opens a new challenge in obtaining nanometric devices with particular optical features suitable for applications. In our case, two different systems have been realized, that are able to work at two or three different wavelengths, depending on the thickness of the thick-ITO layer. The first nanocavity (t1 = 230 nm) is able to confine two different working wavelengths, λ = 390 nm and 550 nm, while the second one (t2 = 350 nm) shows a three wavelengths confinement (λ = 410 nm, 520 nm and 750 nm). We indicated these systems as two- and three-bands metamaterials (2BMM and 3BMM, respectively). The two samples were fabricated by starting with a 20 nm Ag layer deposited on a glass substrate. Then, different thick-ITO slabs (t1 and t2) were deposited. Light confinement was achieved by putting a final Ag/thin-ITO (20/20) bilayer, of 40 nm total thickness, on top of the thick-ITO slabs, obtaining the final 2BMM and 3BMM systems, respectively. The two configurations presented unique optical behaviors that were also verified by the numerical model. The experimental analysis was carried out as a function of the incident angles (50°, 60° and 70°), under which the p-pol reflectances were measured (Fig. 5c and d refer to the 2BMM and the 3BMM, respectively). As the incident angle increases, a blue shift of the reflectance dips is observed, due to the gap surface plasmons (gsp) that are able to establish guiding modes in the ITO thick slabs with t1 (2BMM) and t2 (3BMM) thicknesses. In the presence of an optical cavity, confined modes satisfy the general Fabry–Perot (FP) condition βtcav = mπ − ϕ, where β is the complex wave-vector of the lightwave propagating within the cavity, m identifies the given mode generated inside the cavity, ϕ is the reflection dephase angle and tcav is the cavity thickness. In this case, the FP relation, applied to a nanocavity confined by metal layers,[52] is modified as in the following, βgsptcav = mπ − ϕ, where βgsp is calculated by considering the dispersion relation of a general MIM structure:with , where the subscript m and d refer to the metal and dielectric layers, respectively. The evaluation of βgsp in our cases is reported in the ESI† with the related dispersion graphs (Fig. S2a and c†) and the surface plasmon modes calculations (Fig. S2b and d†).[53-61] In Fig. S2b (ESI†), the modal analysis is reported for the 2BMM, from which it is possible to extract the observed modes. The obtained modes confirm the results shown in Fig. 5a: for m ∼ 3, λ = 390 nm, while λ = 550 nm corresponds to the mode m ∼ 2. For the 3BMM (tcav = 350 nm), the results derived from Fig. S2d† are in agreement with those shown in the electric field maps of Fig. 5b: m ∼ 4 is related to λ = 410 nm, m ∼ 3 to λ = 520 nm and m ∼ 2 to λ = 750 nm. These confined modes are depicted as black curves on the electric field distribution maps reported in Fig. 5a and b. A qualitative confirmation of this fact is given by the spectral position of the dips present in the reflectance curves (Fig. 5c and d) that correspond to the working wavelengths mentioned above, for the 2BMM and the 3BMM, respectively. The curves shown in Fig. 5e and f highlight a value of transmitted light overcoming, for the first system, 60% at λ = 550 nm in the cases of both experiment and simulation. The second cavity instead presents a larger difference at λ = 520 nm between the experimental and simulated systems, even if the transmitted signal still reaches a noticeable 50% value in the real case. We underline that, in specific spectral windows, the presented systems show transmittance values proportionally much higher than in a single Ag layer (20 nm). However, it should be noted that, passing through a nanocavity, light runs into two metal layers (20 nm each) and two (thick and thin) ITO slabs, corresponding to hundreds of nanometers (Fig. 2a, blue squares). The surprisingly high signal transmitted by these asymmetric nanocavities thus evidences their optical behavior as originating from an unusual effective medium. The slight amplitude mismatch between the experimental and numerical reflectance dips and transmittance peaks observed in Fig. 5d and f (3BMM) is probably due to a collapse of the ITO slab in the Ag layer, causing an overall reduction of their effective thicknesses. The field maps and reflectance for s-pol, as well as the comparison of Ψ and Δ between the experiments and the numerical simulations, are reported in the ESI (Fig. S3†); the reflectance curves are also characterized by the same number of dips for this “out-of-plane” polarization. The trend of Ψ presents two or three sigmoidal features corresponding to the number and position of dips in reflectance for both polarizations (Rp and Rs). The same behavior is present in the phase-shift reported in the Δ curves (Fig. S3e and f†).
Fig. 4
Reflectances (Rp, and Rs, red and black lines with symbols) and transmittance (T, blue lines with squares) measured (a) and calculated using the NEA model (b) for the metamaterial with a three-material unit cell (AZO/ITO/Ag). In the inset, a 3D sketch of the sample is reported. Comparison between experimental (colored lines with symbols) and numerical systems (solid lines) Ψ (c) and Δ (d), for the same sample.
Fig. 5
Electric field maps distribution for the 2BMM (a) and 3BMM (b) systems. The black solid lines marked on the field maps represent the even and odd modes supported by the nanometric cavities. Rp reflectances at different angles of incident light, measured and calculated using the NEA model for the 2BMM and the 3BMM, are reported in (c) and (d) respectively. Comparison between experimental and numerical transmittances for the 2BMM and the 3BMM are reported in (e) and (f) respectively.
It is worth noting that the interesting simultaneous presence of dips in the p- and s-pol reflectance curves shown by this MIM structure are independent from the impinging light polarization direction. As such, these results can be exploited to realize sensors or photovoltaic cells that are able to work in a wide UV-VIS-NIR range.[48,62-67] In order to complete the numerical analysis of the 3BMM structure, we performed a parametric study of reflectance and transmittance of the system as a function of wavelength (λi) and incidence angle (θi) of the incoming EM wave in the visible range. The resulting reflectance map, shown in Fig. 6a, reveals the coupling between the impinging light and the 3BMM, evidenced by the three permitted wavelength bands for incident angles ranging between 0° and 85°. Position and number of these bands are strongly related to the thickness of the cavity. For instance, if the thickness is reduced to about 200–300 nm, the dip located around λ = 750 nm disappears, as shown by the 2BMM (Fig. 5c); for a thickness of about 400 nm the wavelength of the same dip is instead shifted to around λ = 800 nm. Fig. S4 of the ESI† demonstrates how sensitive the spectral position is of a given reflectance dip of the 3BMM as a function of the cavity thickness. By considering an incident angle of θi = 50° and varying the cavity thickness from 218 nm to 242 nm, it is possible to move the reflectance dip by Δλ = 40 nm, from 535 nm to 575 nm. Fig. S4b† reports the linear trend of the cavity working wavelength as a function of the cavity thickness. By keeping the cavity thickness fixed, the spectral position of the reflectance dips can also be influenced by the presence of a material on top of the 3BMM surface. This suggests the possible utilization of this system as a sensor for refractive index variations. In order to verify this possibility and carry out a further performance check on the developed tool, we numerically designed another 3BMM structure with 10 nm of AZO placed on its top surface. For this test, the more sensitive confined mode (m ∼ 2 at λ = 750 nm) was selected. By impinging at the Brewster angle (60°, white dashed line in Fig. 6a) and with p-polarized light on the structure, the reflectance curve obtained in the presence of the AZO coating (dashed-dot red line in Fig. 6b) shows a clearly evident redshift of about 10 nm with respect to the curve for the system without the coating (black solid line in Fig. 6b). As a final remark, the exploitation of a generic N-BMM structure for applications presents many benefits: (i) the position and amplitude (depth) of the transmittance (reflectance) are directly correlated to the thickness of the ITO slab; (ii) these transmittance peaks (reflectance dips) have been obtained without the presence of 3D structures on the material surface to couple the impinging light in; (iii) the system is independent from the impinging light polarization; (iv) the fabrication process used to realize the structures is fast and cost-effective.
Fig. 6
(a) Reflectance map calculated by varying the incident wavelength λi and angle θi. The dashed white line drawn on the field map indicates the position of the Brewster angle for this structure that has been used as the incident angle in the numerical simulations whose results are reported in the graph. (b) Rp reflectances calculated for the sensor realized by using the 3BMM with and without an AZO coating layer of 10 nm on its top surface (solid black and dashed-dot red lines, respectively).
Conclusions
In this contribution, we have presented a comprehensive optical analysis of different nanoscale structures. We propose a solid numerical tool, based on the COMSOL Multiphysics platform, that is able to predict the optical behavior of real case studies. As a crucial test-bed of the tool, optical features such as reflectance, and transmittance, as well as the ellipsometric angles Ψ and Δ have been calculated in the presence of systems with increasing structural complexities, from single metal layers to asymmetric multi-layers behaving as optical nanocavities for light. The latter configuration shows a peculiar metamaterial behaviour originating from its unusual effective medium properties, that can be efficiently exploited in sensors applications. Indeed, in specific spectral windows, the asymmetric multi-layer with a full thickness of hundreds of nanometers shows a hyper transmission of light with respect to a single 20 nm-thick Ag layer. The resulting comprehensive analysis of the considered nanoscale systems has shown an excellent agreement between numerical and experimental curves, as detailed by the comparative tables reported in the ESI.† This confirms the effectiveness of the tool as a significant instrument for nanotechnology design and fundamental research in nanoscience.
Authors: Kandammathe Valiyaveedu Sreekanth; Yunus Alapan; Mohamed ElKabbash; Efe Ilker; Michael Hinczewski; Umut A Gurkan; Antonio De Luca; Giuseppe Strangi Journal: Nat Mater Date: 2016-03-28 Impact factor: 43.841
Authors: Antonio Ferraro; Mauro Daniel Luigi Bruno; Giuseppe Papuzzo; Rosa Varchera; Agostino Forestiero; Maria Penolope De Santo; Roberto Caputo; Riccardo Cristofaro Barberi Journal: Nanomaterials (Basel) Date: 2022-04-09 Impact factor: 5.719