Numerical investigation of a refractive index SPR D-type optical fiber sensor using COMSOL multiphysics

Recently, many programs have been developed for simulation or analysis of the different parameters of light propagation in optical fibers, either for sensing or for communication purposes. In this paper, it is shown the COMSOL Multiphysics as a fairly robust and simple program, due to the existence of a graphical environment, to perform simulations with good accuracy. Results are compared with other simulation analysis, focusing on the surface plasmon resonance (SPR) phenomena for refractive index sensing in a D-type optical fiber, where the characteristics of the material layers, in terms of the type and thickness, and the residual fiber cladding thickness are optimized.


Introduction
To estimate the behavior of an optical fiber sensor, it is very important to use a simulation tool to analyze parameters such as: magnetic and electric field intensities, effective refractive index, among others. Several difficulties exist in developing a good simulation program, including the necessary approximations when writing the code for 2D or 3D [1]. In particular, when the sensor structure is complex, the calculation becomes too cumbersome, and it is necessary to use simplified methods, for example: expansion and propagation method (MEP) and the method for multilayer structure transfer matrix modeling [2][3]. The latter allows a better approximation to the optical fiber cylindrical structure [4], but fails to get good results for nano-structures in optical fibers [1]. One solution is to use finite-difference time domain (FDTD), which allows computing the magnetic and electric field distribution, but requires huge quantity of computing memory [5]. In this paper, we demonstrate another method of studying the behavior of optical fiber sensors based on the surface plasmon resonance (SPR) using COMSOL Multiphysics, a commercial program that uses finite element method (FEM).
On the other hand, D-type fiber is an optical fiber with numerous applications in optical sensing for different areas of engineering. Gas detection [6] and curvature sensing [7] are examples of sensing applications using this fiber. The D-type fiber has also been implemented as a biosensor using the surface-plasmon resonance technology [8,9]. In this simulation work, we improved the design of an SPR refractive index sensor, based on a D-type optical fiber, where the characteristics of the material layers, in terms of type and thickness, and the residual fiber cladding thickness were optimized.

Theory
We have analyzed an optical fiber sensor based on SPR composed by a D-type fiber spliced between two single mode fibers, as shown in Fig. 1. The main design parameters of the sensor included the length of the sensor L, the radius and refractive index of the core (r c and n c , respectively) and of the cladding (r cladd and n cladd , respectively), the distance between the core and the metal d (the residual cladding), the thickness of the metal d m , and the refractive index of the external medium n ext . The losses of the light propagating in the fiber were determined by the tuning between the wavelength of the light beam and the SPR which was strongly dependent on the refractive index of the external medium. The transmission coefficient of the sensor could be used to assess with accuracy the value of n ext .

Calculated transmission coefficient using COMSOL
In this work, we have conducted a 2D analysis of the mode structure and the electromagnetic field modes along the transverse plane of a D-type fiber using the mode analysis utilities of COMSOL Multiphysics. The electromagnetic fields in optical fiber waveguides are governed by the macroscopic Maxwell's equations in the absence of currents or external electric charges: where E, H, D and B are the electric, the magnetic, the dielectric and the magnetic induction fields, respectively. Also, the term J is the current density,  is the charge density, r is the spatial coordinate, and t denotes time. The time-harmonic solutions describing strictly monochromatic fields are of the form where  is the angular frequency of light. In this representation, the fields are complex quantities whose real parts correspond to the physical fields [10].
In linear, isotropic, and nonmagnetic media, the following constitutive relations are where 0  and 0  are the permittivity and permeability of free space, respectively, and  r denotes the relative, material-dependent permittivity. In general, these quantities are functions of the spatial coordinates. Taking the curl of (1) and using (3), (7) and (8) the same is the magnetic field: In optical fibers, the dependency in the spatial z coordinate along the axis is obtained using the variable separation method: where i  is the propagation constant of the i-th mode, and r  is the position vector in the plane perpendicular to the optical axis. The solutions of (9) and (10) were obtained using the FEM, which basically consists in dividing the simulation domain into smaller subdomains forming a mesh as shown in Fig. 2. The subdomains have different sizes and are smaller near the interfaces between different media to account for steeper variations of the field. The field equations are then discretized into an algebraic system of equations and solved for their characteristic eigenvalues. The electric and magnetic fields are dependent on the angular frequency [ (9) and (10)] as well as on the refractive index of the materials. For that, it is necessary to calculate the material's refractive index for all frequencies under study. For a dielectric layer, in this case of an optical fiber, it is possible to use the Sellmeier equation [3] where B q , C q are the Sellmeier coefficients, determined experimentally [11], for a germanium-doped silica core fiber and fluorine-doped silica cladding. For the metallic layer, the permittivity and the refractive index can be obtained from the Drude model as where n m and k m are the real and imaginary values of the refractive index for the metal, respectively, m  is the permittivity complex metal, and λ c and λ p denote the plasma wavelength and the collision wavelength, respectively, that is defined in [3] for the metals Ag, Au, Cu and Al. Figure 3 illustrates the intensity of the electric field using the study "Mode Analyses" from COMSOL Multiphysics, for the structure of Fig. 1 and the mesh of Fig. 2 with [Fig. 3(a)] and without [ Fig. 3(b)] a metallic layer, in this case of a 65-nm-thickness gold layer. The wavelength of the light was 802 nm. The study allowed confirming the single-mode behavior propagation in the optical fiber. The 1D electrical field amplitude in the optical fiber is also shown in Figs. 4(a) and 4(b), with and without a metallic layer, respectively. Comparing both figures, it is possible to see the electric field intensity external to the fiber is stronger when using the metal.

Transmission coefficient calculated using Fresnel laws
To verify that the method works properly, a comparison was made with the implemented algorithm in [12], which used Fresnel equations applied to the structure in Fig. 1 where k i is the component of the wave vector of the interface of the two layers of the sensor in the direction z and is given as where n 1 , n 2 , n 3 and n 4 represent the refractive indices of the core, cladding, metal and the external test medium, respectively.
Applying (16) and (17), it is possible to obtain the results by the two different methods, as shown in Fig. 5. The behavior of the two methods was similar, having a difference between the transmission coefficients and a small shift in the wavelength dips. We attributed this difference to the fact that in the Fresnel equations' algorithm it is only considered planar waves in a fairly symmetrical arrangement. On the other hand, when using the FEM, we considered the D-type fiber as a non-symmetrical cylindrical waveguide, being able to model the inhomogeneous optical regions with a resolution of the cell size, resulting in a more accurate outcome. In terms of the results and in what concerns the material thickness, the optimal point occurred when a layer with a thickness of 55 nm to 65 nm was used.
Another way to test the efficiency of the described procedure is to calculate the sensor sensitivity (/RIU) as function of the different refractive indices of the external environment. The sensitivity for both methods is almost equal and is close to 3150 nm/RIU [2].
Tailoring the simulation analysis in COMSOL Multiphysics, it is possible to optimize the sensitivity, transmission coefficient dip, wavelength operation area, amongst others for a refractive index SPR D-type optical fiber sensor. To decrease the depth of the transmission coefficient dip and consequently lower the sensitivity of the external medium, d can be increased, as shown in Fig. 6. From Figs. 5 and 6, it is possible to have a sensor that works for an area of operation near 820 nm, for a metal thickness of 65 nm (Au).
In case another wavelength is required, one possible solution is to apply an additional layer of a dielectric with a high refractive index, such as tantalum pentoxide (Ta 2 O 5 ) [1][2], which simulation results can be seen in Fig. 7 and compared with the results presented in [2]. For different thicknesses of Ta 2 O 5 , the transmission coefficient dip of the sensor operation is not significantly altered being possible to tailor the wavelength sensor operation [2].

Conclusions
Another method was demonstrated to study the behavior of optical fiber sensors for refractive index measurement based on SPR using COMSOL Multiphysics, a commercial program that uses the finite element method. The two simulations, one with COMSOL Multiphysics and the other with Fresnel's equations, present a similar behavior. The graphical interface of COMSOL Multiphysics facilitates the simulation work, having no need to develop complex formulations. Also, it has the ability to model inhomogeneous optical regions with a resolution of the cell size and allows the analysis of other parameters such as the intensity of magnetic and electric field across the structure [13]. COMSOL Multiphysics permits in a graphical environment more accurate and realistic results than traditional approaches, although at the expense of longer running time.
It was also possible to demonstrate the use of COMSOL Multiphysics to improve the performance of a refractive index SPR D-type optical fiber sensor, where the characteristics of the material layers, in terms of the type and thickness, and the residual fiber cladding thickness are optimized.
Open Access This article is distributed under the terms of the Creative Commons Attribution License which permits any use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.