Second-Order Dispersion Sensor Based on Multi-Plasmonic Surface Resonances in D-Shaped Photonic Crystal Fibers

: This paper proposes a scheme to determine the optical dispersion properties of a medium using multiple localized surface plasmon resonances (SPR) in a D-shaped photonic crystal ﬁber (PCF) whose ﬂat surface is covered by three adjacent gold layers of different thicknesses. Using computational simulations, we show how to customize plasmon resonances at different wavelengths, thus allowing for obtaining the second-order dispersion. The central aspect of this sensing conﬁguration is to balance miniaturization with low coupling between the different localized plasmon modes in adjacent metallic nanostructures. The determination of the optical dispersion over a large spectral range provides information on the concentration of different constituents of a medium, which is of paramount importance when monitoring media with time-varying concentrations, such as ﬂuidic media.


Introduction
Surface plasmon polaritons are electromagnetic modes that arise from the coupling between photons and free-electron oscillations at a conducting surface [1]. The practical excitation of a surface plasmon resonance (SPR) promotes the confinement of optical power at subwavelength dimensions and represents a milestone in the development of nanooptical sensors. Since surface plasmon modes are quite dependent on the refractive index of the surrounding medium, higher levels of sensitivity and resolution in a broad spectral range can be reached, which is the central reason for the huge effort dispended in related studies in the two last decades [2]. As a result, the losses variations and mode phase shifts can be used to describe medium properties in terms of its refractive index [3].
A conventional SPR prism-based sensor, although rather efficient to excite surface plasmon modes [4,5], suffers drawbacks due to its bulky size, which hampers remote sensing applications, and high-cost fabrication process. As an alternative, a common configuration consists of an optical fiber with partial cladding removed for the deposition of a thin metallic layer. The plasmonic resonance is achieved when the fiber core guided mode phase is equal to the surface plasmon mode phase at the conducting surface. Multiple sensing devices are designed to reach this condition as tapered fibers, grating-based fiber sensors and D-shaped fibers. In [6], a tapered optical fiber covered by gold nanoparticles is applied for biomolecular sensing. By its turn, the deposition of gold nanoparticles over a grating fiber sensor provides a highly sensitive and selective platform, which is quite useful for chemical applications [7,8]. In [9], a surface plasmon resonance sensor operating at the near-infrared band is constructed by coating the flat face of a D-shaped optical fiber with a graphene-based metal oxide layer.
Photonic crystal fibers (PCF) are highly suitable to excite and enhance surface plasmon resonances due to their unique characteristics, such as fine control of the evanescent field penetration into the conducting medium and high mode confinement with a large mode area [10][11][12]. A PCF has an arrangement of air holes periodically distributed over its crosssectional area. The light waveguiding is favored by inserting a defect in the set of air holes. The modal properties of the fiber are extremely influenced by the structural parameters, such as the diameter of the air holes and the distance between them. Several designs of SPR sensors based on D-shaped PCF with high levels of sensitivity and resolution have been reported [12][13][14][15]. However, to the best of our knowledge, all these sensors are used to interrogate the average refractive index or a specific wavelength of the surrounding medium, with no further information regarding the wavelength dependence of the analyte.
We reported a D-shaped PCF with two gold slabs [16] whose response has a highly linear dependence on the plasmonic resonances with the first-order dispersion. The designed dispersion sensor finds potential applications for the investigation of changes in fluidic media over time, for instance. This is because small variations on the chemical composition produce subtle changes in the dispersion relation of the medium [17,18], which cannot be roughly monitored by conventional refractive-index-based sensors.
In this work, we investigate theoretically the viability of applying an SPR sensor based on a D-shaped PCF to monitor the second-order dispersion of an optical medium. The design is projected to provide three distinct and independent plasmonic resonances. This is obtained by depositing three gold slabs with different thicknesses on the top of the core region at the flat surface of the PCF. The sensing response is taken by the changes in the amplitude of the resonance peaks and the distance between them for parabolic dispersive profiles. Our final goal is to demonstrate a relationship between the sensor response with the input second-order dispersive medium. In the following sections, we detail the modeling of the SPR sensor based on a D-shaped PCF and discuss the results for second-order dispersion sensing.

Designed Structure and Modeling
The cross-section of the proposed sensor is shown in Figure 1. It is composed of a hexagonal arrangement of air holes in silica background. The fiber core (the central region with an absence of air holes) has a diameter of 5 µm, the pitch (distance between two adjacent air holes) is Λ = 2 µm and the diameter of each air hole is related to the pitch by d/0.88 = Λ. This type of D-shaped structure can be obtained by the stack-and-draw process [19] followed by a side-polishing or controlled etching technique [20]. Three gold slabs of equal width w = 1.5 µm but with different thicknesses (t 1 = 40 nm, t 2 = 15 nm and t 3 = 30 nm) are deposited on the flat surface of the fiber, which can be performed in practice by the CVD (Chemical Vapor Deposition) process [21].
The Finite-Element-Method-based software COMSOL Multiphysics [22] is applied for the numerical modeling. The computational domain comprises the cross-section of Figure 1a with diameter D = 24 µm, truncated by a 0.1D thick PML (Perfectly Matched Layer). More specifically, the Wave Optics package in the frequency domain is applied to carry out 2D computational simulations to obtain the eigenvalues of the Helmholtz equation in the angular frequency ω: where k 0 is the magnitude of the free-space wavenumber and ε is the complex frequencydependent relative permittivity. E(r ⊥ , ω) is the modal electric field distribution at the position r ⊥ perpendicular to the direction of light propagation. The expression for ε depends on the region at which (1) is solved. The complex effective index n e f f is the result of the modal fields that arise from the coupling between the fundamental fiber mode and the plasmonic excitations at the boundaries of the gold layers. While the real part of n e f f refers to the mode phase, the imaginary part is related to the losses experienced by the confined mode as the fields are tunneled through the gold slabs. The Finite-Element-Method-based software COMSOL Multiphysics [22] is applied for the numerical modeling. The computational domain comprises the cross-section of Figure 1a with diameter D = 24 µm, truncated by a 0.1D thick PML (Perfectly Matched Layer). More specifically, the Wave Optics package in the frequency domain is applied to carry out 2D computational simulations to obtain the eigenvalues of the Helmholtz equation in the angular frequency : where is the magnitude of the free-space wavenumber and is the complex frequency-dependent relative permittivity.
, is the modal electric field distribution at the position perpendicular to the direction of light propagation. The expression for depends on the region at which (1) is solved. The complex effective index is the result of the modal fields that arise from the coupling between the fundamental fiber mode and the plasmonic excitations at the boundaries of the gold layers. While the real part of refers to the mode phase, the imaginary part is related to the losses experienced by the confined mode as the fields are tunneled through the gold slabs.
For accurate modeling, the dispersive character of all materials in the sensing structure has to be considered. The refractive index of silica is computed from the Sellmeier equation [23] and the air region is modelled by a constant refractive index nair = 1. Drude-Lorentz formalism [24], which is an improvement on the original Drude model, is used to characterize the dispersive character of gold. This model provides a better approach to the corresponding experimental data: The parameters ωp = 2155.6 THz and γ = 15.92 THz are the plasma frequency and damping factor, respectively, and ε∞ = 5.9673 stands for the residual polarization of gold at high frequencies. The parameters in the correction factor are Δϵ = 1.09, ΩL = 650.07 THz and ΓL = 104.86 THz [24].
The modeling approach allowed for solving (1) without the requirements of significative computational efforts and extensive time consumption. For accurate modeling, the dispersive character of all materials in the sensing structure has to be considered. The refractive index of silica is computed from the Sellmeier equation [23] and the air region is modelled by a constant refractive index n air = 1. Drude-Lorentz formalism [24], which is an improvement on the original Drude model, is used to characterize the dispersive character of gold. This model provides a better approach to the corresponding experimental data: The parameters ω p = 2155.6 THz and γ = 15.92 THz are the plasma frequency and damping factor, respectively, and ε ∞ = 5.9673 stands for the residual polarization of gold at high frequencies. The parameters in the correction factor are ∆ = 1.09, Ω L = 650.07 THz and Γ L = 104.86 THz [24].
The modeling approach allowed for solving (1) without the requirements of significative computational efforts and extensive time consumption.

Sensitivity Performance Analysis
In the sensing platform schematically depicted in Figure 1, the resonance condition is achieved when there is phase-matching between the fundamental fiber mode and surface plasmon mode. Due to this phenomenon, the spectral confinement loss presents a peak, which is related to the fact that the amount of modal energy penetrating into the gold layers is maximal. When the refractive index of the analyte changes, the wavelength that satisfies the resonance condition also changes. In Figure 2, the imaginary part, expressed in terms of the spectral losses, and the real part of the effective indexes of the fundamental fiber mode and plasmonic modes are shown for an analyte with dispersive parabolic refractive index Ri analyte = 7.41 × 10 −6 λ 2 − 0.0083λ + 3.65, with the wavelength λ measured in nanometers. The effective indexes were computed from the numerical modeling described previously applied to the solution of (1) over the entire sensing structure. At the resonance wavelengths, the effective indexes of the surface plasmon mode and fiber mode present the same real part due to phase-matching and, as a consequence, the spectral losses of the fiber mode reach a peak. These modes are uncoupled at non-resonance wavelengths, and the corresponding effective indexes are completely distinct. By adding metallic layers of different thicknesses, several resonant peaks are expected in the spectral losses. In Figure 2, the real part of the effective index of the fundamental mode intercepts SPP modes at the wavelengths 568, 586 and 632 nm, respectively. Even though the fundamental mode of the PCF is degenerate for two orthogonal polarizations, only the polarization perpendicular to the gold interfaces allows for the excitation of plasmonic modes.
peak, which is related to the fact that the amount of modal energy penetrating into the gold layers is maximal. When the refractive index of the analyte changes, the wavelength that satisfies the resonance condition also changes. In Figure 2, the imaginary part, expressed in terms of the spectral losses, and the real part of the effective indexes of the fundamental fiber mode and plasmonic modes are shown for an analyte with dispersive parabolic refractive index = 7.41 × 10 − 0.0083 + 3.65, with the wavelength λ measured in nanometers. The effective indexes were computed from the numerical modeling described previously applied to the solution of (1) over the entire sensing structure. At the resonance wavelengths, the effective indexes of the surface plasmon mode and fiber mode present the same real part due to phase-matching and, as a consequence, the spectral losses of the fiber mode reach a peak. These modes are uncoupled at non-resonance wavelengths, and the corresponding effective indexes are completely distinct. By adding metallic layers of different thicknesses, several resonant peaks are expected in the spectral losses. In Figure 2, the real part of the effective index of the fundamental mode intercepts SPP modes at the wavelengths 568, 586 and 632 nm, respectively. Even though the fundamental mode of the PCF is degenerate for two orthogonal polarizations, only the polarization perpendicular to the gold interfaces allows for the excitation of plasmonic modes. In [16], we have explored the dependence of SPR spectra on the first-order dispersion profile of the refractive index of the analyte. In the present work, we focus on the relation of the second-order dispersion with the spectral response of the sensor. For that purpose, we carried out computations of the SPR spectral losses for three dispersive characters: constant refractive index, normal dispersion and anomalous dispersion. The spectral losses of the SPR D-shaped PCF with three gold slabs and their dispersion regimes are shown in Figure 3. For a constant refractive index equal to 1.36 (solid red line), three distinct regions of maximum absorption appear (solid blue line) due to the nonuniformity in the thickness of the three gold slabs at the flat face of the sensor.
Moreover, we consider two arbitrary media with parabolic dispersion profiles of opposite concavities. In Figure 3, the dotted blue line represents the spectral losses due to the anomalous dispersion curve (dotted red line), while the dashed blue line is related to In [16], we have explored the dependence of SPR spectra on the first-order dispersion profile of the refractive index of the analyte. In the present work, we focus on the relation of the second-order dispersion with the spectral response of the sensor. For that purpose, we carried out computations of the SPR spectral losses for three dispersive characters: constant refractive index, normal dispersion and anomalous dispersion. The spectral losses of the SPR D-shaped PCF with three gold slabs and their dispersion regimes are shown in Figure 3. For a constant refractive index equal to 1.36 (solid red line), three distinct regions of maximum absorption appear (solid blue line) due to the nonuniformity in the thickness of the three gold slabs at the flat face of the sensor.
Moreover, we consider two arbitrary media with parabolic dispersion profiles of opposite concavities. In Figure 3, the dotted blue line represents the spectral losses due to the anomalous dispersion curve (dotted red line), while the dashed blue line is related to the normal dispersion curve (dashed red line). Just as in the case of the constant refractive index, there are three distinct resonance wavelengths for a dispersive index. It is also observed that the magnitude of the losses and the distance between resonance wavelengths is influenced by the dispersion profile.
In a comparative analysis in the case of a constant refractive index, the distances between the wavelengths of the three peaks decrease for the normal dispersion and the amplitude losses increase at 582 and 593 nm. On the other hand, the anomalous dispersion curve leads to a larger difference between resonance wavelengths but to smaller amplitude losses at 568 and 586 nm. For the resonance peak at 632 nm, losses remain close for the different dispersion profiles. In a comparative analysis in the case of a constant refractive index, the distances between the wavelengths of the three peaks decrease for the normal dispersion and the amplitude losses increase at 582 and 593 nm. On the other hand, the anomalous dispersion curve leads to a larger difference between resonance wavelengths but to smaller amplitude losses at 568 and 586 nm. For the resonance peak at 632 nm, losses remain close for the different dispersion profiles.

Sensing Response at Second-Order Dispersion
To access the ability of the sensor in retrieving the second-order dispersion of a medium, we applied the multivariate polynomial regression model [25] to obtain the coefficients a, b and c of the general quadratic expression = + + from nine parameters obtained from computed SPR spectra, namely: the amplitudes of the three SPR peaks, the corresponding resonance wavelengths and the distances between these wavelengths. The numerical computation of the SPR spectrum was carried out for 270 distinct parabolic dispersion curves in order to achieve a dataset able to ensure an accurate fitting. The multivariate polynomial regression equation to achieve the coefficients βn can be written in the form as [25]: where y is the vector containing the parabolic profile parameters (specifically, the a's, b's and c's), are the multivariate polynomial regression coefficients, ′ are the vectors containing the simulated parameters, specifically, the resonance amplitude for the three peaks, the wavelength peak for the three resonances and the product among themselves. The independent term ε is related to the errors.
As an example, Figure 4 shows the points retrieved by applying the multivariate polynomial regression to the SPR spectral parameters for both dispersion regimes.

Sensing Response at Second-Order Dispersion
To access the ability of the sensor in retrieving the second-order dispersion of a medium, we applied the multivariate polynomial regression model [25] to obtain the coefficients a, b and c of the general quadratic expression Ri analyte = aλ 2 + bλ + c from nine parameters obtained from computed SPR spectra, namely: the amplitudes of the three SPR peaks, the corresponding resonance wavelengths and the distances between these wavelengths. The numerical computation of the SPR spectrum was carried out for 270 distinct parabolic dispersion curves in order to achieve a dataset able to ensure an accurate fitting. The multivariate polynomial regression equation to achieve the coefficients β n can be written in the form as [25]: where y is the vector containing the parabolic profile parameters (specifically, the a's, b's and c's), β n are the multivariate polynomial regression coefficients, x n s are the vectors containing the simulated parameters, specifically, the resonance amplitude for the three peaks, the wavelength peak for the three resonances and the product among themselves. The independent term ε is related to the errors.
As an example, Figure 4 shows the points retrieved by applying the multivariate polynomial regression to the SPR spectral parameters for both dispersion regimes.
To characterize quantitatively the performance of the fitting model, we computed the error in estimating the parabolic coefficients a, b and c (Ri analyte = aλ 2 + bλ + c) from the multivariate linear regression model. The histograms depicted in Figure 5 show the Gaussian distribution of the error, which is centered at 0 for the three coefficients. To compare the quality of the regression model for future works, we also calculated the fullwidth half-maximum for the histograms: 0.1503, 0.0894 and 0.1521 for the parameters a, b and c, respectively. Moreover, the calculated determination coefficients R 2 = 0.908129 for a, R 2 = 0.905276 for b and R 2 = 0.901035 for c demonstrate a reasonable level of accuracy for the estimative of the parabolic dispersion curve from the SPR spectral parameters for the proposed sensing platform.  To characterize quantitatively the performance of the fitting model, we computed the error in estimating the parabolic coefficients a, b and c ( = a + + ) from the multivariate linear regression model. The histograms depicted in Figure 5 show the Gaussian distribution of the error, which is centered at 0 for the three coefficients. To compare the quality of the regression model for future works, we also calculated the full-width half-maximum for the histograms: 0.1503, 0.0894 and 0.1521 for the parameters a, b and c, respectively. Moreover, the calculated determination coefficients R 2 = 0.908129 for a, R 2 = 0.905276 for b and R 2 = 0.901035 for c demonstrate a reasonable level of accuracy for the estimative of the parabolic dispersion curve from the SPR spectral parameters for the proposed sensing platform.
(a)  To characterize quantitatively the performance of the fitting model, we computed the error in estimating the parabolic coefficients a, b and c ( = a + + ) from the multivariate linear regression model. The histograms depicted in Figure 5 show the Gaussian distribution of the error, which is centered at 0 for the three coefficients. To compare the quality of the regression model for future works, we also calculated the full-width half-maximum for the histograms: 0.1503, 0.0894 and 0.1521 for the parameters a, b and c, respectively. Moreover, the calculated determination coefficients R 2 = 0.908129 for a, R 2 = 0.905276 for b and R 2 = 0.901035 for c demonstrate a reasonable level of accuracy for the estimative of the parabolic dispersion curve from the SPR spectral parameters for the proposed sensing platform.

Conclusions
The characterization of the optical dispersion of the refraction index over a large spectral range can be useful in the determination of the concentration of different constituents of a medium, which is particularly relevant when monitoring media in real time. By adding three gold slabs of different thicknesses at the flat face of the D-shaped photonic crystal fiber, we were able to characterize the second-order dispersion of a medium from the excitation of multiple surface plasmon resonances. This new approach is able to determine the non-local spectral character of the refraction index of realistic media, such as fluids whose composition varies with time. Moreover, since the surface plasmon resonances are excited at distinct spectral channels, the sensing structure can be used to determine simultaneously more than one parameter. Funding: This study was financed in part by the Coordenação de Aperfeiçoamento de Pessoal de Nível Superior-Brasil (CAPES)-Finance code 001. During the period of elaboration of this work, the authors have obtained support from the Conselho Nacional de Desenvolvimento Cientifico e Tecnologico (CNPq). It was also financed by the ERDF-European Regional Development Fund

Conclusions
The characterization of the optical dispersion of the refraction index over a large spectral range can be useful in the determination of the concentration of different constituents of a medium, which is particularly relevant when monitoring media in real time. By adding three gold slabs of different thicknesses at the flat face of the D-shaped photonic crystal fiber, we were able to characterize the second-order dispersion of a medium from the excitation of multiple surface plasmon resonances. This new approach is able to determine the non-local spectral character of the refraction index of realistic media, such as fluids whose composition varies with time. Moreover, since the surface plasmon resonances are excited at distinct spectral channels, the sensing structure can be used to determine simultaneously more than one parameter.