Silicon-Based Multilayer Waveguides for Integrated Photonic Devices from the Near to Mid Infrared

: Advancements in spectroscopy, quantum optics, communication, and sensing require new classes of integrated photonic devices to host a wide range of non-linear optical processes involving wavelengths from the visible to the infrared. In this framework, waveguide (WG) structures designed with innovative geometry and materials can play a key role. We report both ﬁnite element modeling and experimental characterization of silicon nitride multilayer WGs from the visible to the mid-infrared spectral regions. The simulations evaluated optical behavior and mechanical stress as a function of number of WG layers and photonic structure dimensions. WGs were optimized for waveguiding at 1550 nm and 2640 nm. Experimental characterization focused on optical behavior and coupling losses from 532 nm to 2640 nm. Measured losses in WGs indicate a quasi-perfect waveguiding behavior in the IR range (with losses below 6 dB), with a relevant increase (up to 20 dB) in the visible range.


Introduction
In the last few years, Quantum Technologies have been recognized as key for novel disruptive instruments and future industrial products, with the prospect to become pervasive in all application areas, as well as in everyday life [1][2][3]. This "second quantum revolution" will need miniaturized photonic platforms to be embedded in newly designed instruments and products. Therefore, novel laser configurations and photonic components will underpin such developments [4,5] in every spectral region, from the UV to the far infrared [6,7]. In this context, mid-infrared (MIR) can play a key role. Indeed, this is the so-called fingerprint region, already crucial for a number of civilian as well as defense and security applications: from spectroscopy [8][9][10] to free-space communication [11,12], from metrology to imaging [13][14][15][16], and from photonics [17][18][19][20][21] to quantum technologies [22]. In particular, waveguide WG devices are very important for photonics, communication, and health applications [7,[20][21][22][23][24]. MIR WG structures come in a number of geometric shapes and materials, e.g., GaAs/AlGaAs, III-V/silicon, silicon-on-insulator (SOI), silicon-onsapphire, silicon, silicon-on-nitride, silicon-on-calcium-fluoride, silicon-on-lithium-niobate, germanium-on-silicon, or germanium-on-silicon-nitride platforms [18,[25][26][27][28]. However, MIR optical components (detectors, fibers, waveplates, polarizers, optical isolators, etc.) are generally less efficient and have higher production costs than their counterparts in the near infrared (NIR). In this context, it would be desirable to have very broadband WGs working from the MIR to the NIR/visible region. In order to have an integrated platform with components working from visible to MIR, the use of non-linear frequency conversion processes based on silicon is crucial, and exploits the advantages of a well-established technology, providing scalability and cost reduction. In particular, the latest studies report non-linear frequency up-conversion in quasi-phase-matched (QPM) Si-based WGs [27][28][29][30][31][32][33][34][35], particularly second-harmonic generation (SHG) and sum-frequency generation (SFG) [36,37]. Recent experiments with Si-based devices reported non-linear generation with output powers above the nW range [38][39][40]. In other recent studies, silicon nitride (Si 3 N 4 ) was used for SHG [41], although with an efficiency of about two orders of magnitude lower than in traditional non-linear crystals. It is well known that Si 3 N 4 is a centrosymmetric material, therefore lacking any bulk second-order non-linearity; it is nevertheless observed in Si 3 N 4 or Si x N y thin films [42][43][44]. Another well-known possibility with Si 3 N 4 is to strain silicon WGs, thus enhancing the non-linear response [31,45,46], but in this case the non-linearity is mostly due to the silicon itself [47,48]. To date, an extremely small χ (2) coefficient in integrated platforms is attributed to the presence of an asymmetric dipole potential at the interface between the two materials (i.e., Si 3 N 4 and SiO 2 ). In this way, it is possible to break the center-symmetry of the two materials in a limited spatial region [41]. In previous studies, large propagation losses were reported [49]. This was explained by poor optical confinement of guided modes, due to a WG core thickness that was too small. In particular, the large strain occurring between the core and the cladding layers limits the thickness value [40,50].
One of the targets of this study is the evaluation of the optical properties of waveguiding structures in different spectral regions. This is the first step to enabling non-linear processes. Indeed, to host these processes, a platform is needed that would be capable of guiding a wider range of wavelengths. The following steps are a WG design optimizing phase matching and using opportune materials to enhance the desired non-linear processes.
In the present study, we improve optical confinement by testing non-stoichiometric Si x N y /SiO 2 multilayer WGs (hereafter, we refer to the layer material as SiN). These devices were fabricated and provided by the Bruno Kessler Foundation (FBK) and the University of Trento. The presence of intermediate, thin SiO 2 layers was used to reduce the tensile stress fraction of the overall stress acting on the WG. Thus, a WG core with a total SiN thickness larger than the critical one could be obtained. Our study is divided into two parts: In the first part, we perform a finite element modeling (FEM) of multilayer structure WGs focused on waveguiding behavior in the telecom range to host non-linear optical processes; in the second part, we carry out an initial optical characterization of WGs built according to the theoretical modeling.

Geometry and Simulations
We explored a geometry based on non-stoichiometric silicon nitride (SiN) in combination with silicon dioxide (SiO 2 ) ( Figure 1). This choice is motivated by the need to have both a refractive index larger than that of SiO 2 (2.35 vs. 1.46) [48,51], resulting in smaller waveguide dimensions, and a large energy bandgap (close to 5 eV), resulting in lower material losses [44] and the absence of a two-photon absorption (TPA) in the NIR/MIR. We performed two kinds of simulations for the WGs:

2.
Mechanical stress in the multilayer structure. The implications of these simulations concerned possible local second-order non-linearities of such devices. For this study, a superficial charge was applied at layer interfaces, according to the number of layers. The Young's modulus, Poisson coefficient, and density were fixed to values from the relevant literature [54]. The results are given in the next paragraphs.
Waveguides were fabricated by FBK in order to carry out an initial evaluation of their properties. We investigated the losses experienced along the multilayer WG, measuring the transmitted power value at one WG end.  Figure 2 shows the free-space optical setup used for the characterization of the WGs provided by FBK and the University of Trento. The samples available for measurement had a common thickness size (770 nm), with a length of 6 mm and five different widths (7, 8.3, 8.6, 9.4, and 10 µm-see also the last figure for info about stress for every sample). The goal of this experiment was an initial evaluation of the WG properties (measuring the transmitted power value at one WG end) at different wavelengths and to compare those results with the FEM analysis described in next sections.

Experimental Setup
We used four different laser sources to monitor the waveguiding behavior from the visible to the MIR range: a laser diode emitting 0.9 mW at 532 nm (Thorlabs CPS532), a laser diode with 1.2 mW of power at 635 nm (Thorlabs CPS635), an Er-doped fiber laser with 0.8 W at 1550 nm (IPG Photonics), and a DFB laser with 7 mW at 2640 nm (Nanoplus). Schematic of the free-space optical system assembly for the characterization of the waveguiding behaviour (in text and figure, M is for mirror and C is for collimator). Inset: Optical pictures of behaviour at 635 nm of a 6 mm-length WG with a width of 9 µm. Waveguides were fabricated by FBK in order to carry out an initial evaluation of their properties. We investigated the losses experienced along the multilayer WG, measuring the transmitted power value at one WG end. Figure 2 shows the free-space optical setup used for the characterization of the WGs provided by FBK and the University of Trento. The samples available for measurement had a common thickness size (770 nm), with a length of 6 mm and five different widths (7, 8.3, 8.6, 9.4, and 10 µm-see also the last figure for info about stress for every sample). The goal of this experiment was an initial evaluation of the WG properties (measuring the transmitted power value at one WG end) at different wavelengths and to compare those results with the FEM analysis described in next sections. Waveguides were fabricated by FBK in order to carry out an initial evaluation of their properties. We investigated the losses experienced along the multilayer WG, measuring the transmitted power value at one WG end.  Figure 2 shows the free-space optical setup used for the characterization of the WGs provided by FBK and the University of Trento. The samples available for measurement had a common thickness size (770 nm), with a length of 6 mm and five different widths (7, 8.3, 8.6, 9.4, and 10 µm-see also the last figure for info about stress for every sample). The goal of this experiment was an initial evaluation of the WG properties (measuring the transmitted power value at one WG end) at different wavelengths and to compare those results with the FEM analysis described in next sections.

Experimental Setup
We used four different laser sources to monitor the waveguiding behavior from the visible to the MIR range: a laser diode emitting 0.9 mW at 532 nm (Thorlabs CPS532), a laser diode with 1.2 mW of power at 635 nm (Thorlabs CPS635), an Er-doped fiber laser with 0.8 W at 1550 nm (IPG Photonics), and a DFB laser with 7 mW at 2640 nm (Nanoplus).  We used four different laser sources to monitor the waveguiding behavior from the visible to the MIR range: a laser diode emitting 0.9 mW at 532 nm (Thorlabs CPS532), a Appl. Sci. 2021, 11, 1227 4 of 14 laser diode with 1.2 mW of power at 635 nm (Thorlabs CPS635), an Er-doped fiber laser with 0.8 W at 1550 nm (IPG Photonics), and a DFB laser with 7 mW at 2640 nm (Nanoplus).
Our experimental setup ( Figure 2) enabled a fast switch of the different sources and detectors. Flip mirrors were used to easily switch between visible (M7) and IR (flip mirror) input light in the tested WGs. The remaining propagation optics, except WG-input and WG-output collimators were shared at all tested wavelengths by using gold-coated mirrors. For the collimation optics (C1 and C2 in the figure), a couple of collimation lenses were used, anti-reflective coated between 350 and 700 nm for the visible sources and between 1.6 and 3 µm for the IR sources. Input laser beams at the different wavelengths were optically shaped to match the WG' mode field diameters. In Figure 2 lines of different colors are depicted according to the different spectral regions: visible (green line), NIR (red line), and MIR (brown line). Half-wave plates were placed in the NIR and MIR lasers' input paths in order to control beam polarization. Fine realignment of the input and output optics to optimize waveguiding and detection was performed when switching between laser sources.
At the WG output, the light was sent to three different detectors, depending on the source in use. In addition, to analyze, check, and optimize the waveguiding behavior of the samples, as well as to control the correct IR beam line-up, an additional branch could be selected, by means of the flip mirror, after M5 to send the output light to a charge-coupled device (CCD) camera (an Si sensor for visible light, and an InSb camera for NIR/MIR light).

Optical Modeling
In order to guide a wide range of wavelengths from visible to IR, a SiN multilayer WG structure was designed to guide those and host non-linear optical processes between a TEM-polarized pump and signal modes with NIR-MIR wavelengths. The basic WG structure that was studied is shown in Figure 1. The refractive indexes of SiO 2 , SiN, and air were found in the relevant literature [55,56].
We changed the WG dimensions and number of layers to optimize the waveguiding behavior in a wide emission spectrum range, as well as to host the non-linear optical process involving NIR and MIR wavelengths. The cross-section dimension of the surrounding air was selected according to WG dimensions; the distance between the WG edge and the simulation boundary was set to 2.5 µm. In addition, a 4-µm SiO 2 substrate was modeled. That thickness was chosen considering the SOI wafer used for fabrication (the thickness of the SOI wafer before processing was 1 µm for the device layer, 4 µm for the SiO 2 layer, and 650 µm for the substrate).
When working with interlayers, it is convenient to define a transition boundary condition, because a discontinuity in the tangential electric field is clearly present [57]. Transition boundary conditions are characterized by: , k = ω (ε + (σ(jω)))µ (2) where µ (relative permeability), ε (relative permittivity), and σ (electrical conductivity) are all taken from SiN. E T1,2 are the transversal components of the electric field, and k is the wavenumber. The transition boundary thickness (d) is computed as 1 Å.
A tetrahedral grid with a maximum-element size below λ/20 [58] was used to mesh the computational model. Thereafter, the electric field distribution was calculated assuming Maxwell's equations: where E is the electric field and k 0 is the free-space wavenumber. Resolving this equation, we can obtain the effective refractive index given by: Due to the tolerance and robustness of the fabrication method, and based on previous papers [45], a working range of values was defined for the different WG dimensions (i.e., SiN and SiO 2 thickness). Therefore, WG widths and layer thicknesses were calculated, respectively, within the 0.5-20 µm and 50-500 nm size intervals.
A summary of the simulated model behavior is reported in Figure 3, where the effective refractive index (n eff ) of the guided fundamental single mode is plotted vs. the SiO 2 interlayer thickness. The SiO 2 thickness values below 50 nm were strongly affected by mesh size and computation errors. From the initial results, we observed that the waveguiding effectiveness dramatically increased for a SiO 2 thickness between 50 and 150 nm, when the WG is fabricated with more than three SiO 2 interlayers. Actually, a simulation of the model behavior for different WG widths was carried out. where E is the electric field and k0 is the free-space wavenumber. Resolving this equation, we can obtain the effective refractive index given by: Due to the tolerance and robustness of the fabrication method, and based on previous papers [45], a working range of values was defined for the different WG dimensions (i.e., SiN and SiO2 thickness). Therefore, WG widths and layer thicknesses were calculated, respectively, within the 0.5-20 µm and 50-500 nm size intervals.
A summary of the simulated model behavior is reported in Figure 3, where the effective refractive index (neff) of the guided fundamental single mode is plotted vs. the SiO2 interlayer thickness. The SiO2 thickness values below 50 nm were strongly affected by mesh size and computation errors. From the initial results, we observed that the waveguiding effectiveness dramatically increased for a SiO2 thickness between 50 and 150 nm, when the WG is fabricated with more than three SiO2 interlayers. Actually, a simulation of the model behavior for different WG widths was carried out. In Figure 4, a comparative study for the TM00 and TM01 modes at 1550 nm is plotted for a WG with three SiO2 interlayers. In Figure 4a, a quasi-constant neff value for TM00 with WG widths above 4 µm can be observed. However, neff for TM01 becomes constant for WG widths above 10 µm. On the other hand, mode losses were calculated, as we show in Figure 4b. As we can observe, for a reduced lateral dimension (below 3 µm) of the WG crosssection, light confinement is not possible along the WG (Figure 4a), and, therefore, the mode losses are significantly incremented as a result of the scattering losses increase (Figure 4b). In Figure 4, a comparative study for the TM 00 and TM 01 modes at 1550 nm is plotted for a WG with three SiO 2 interlayers. In Figure 4a, a quasi-constant n eff value for TM 00 with WG widths above 4 µm can be observed. However, n eff for TM 01 becomes constant for WG widths above 10 µm. On the other hand, mode losses were calculated, as we show in Figure 4b. As we can observe, for a reduced lateral dimension (below 3 µm) of the WG cross-section, light confinement is not possible along the WG (Figure 4a), and, therefore, the mode losses are significantly incremented as a result of the scattering losses increase (Figure 4b).
Appl. Sci. 2021, 11, x FOR PEER REVIEW 6 of 14 In this context, we have simulated the cut-off wavelength for different interlayer numbers, as shown in Figure 5a. Furthermore, in Figure 5b, the mode loss and the effective index at cut-off wavelengths are plotted. A dramatic increase in mode loss is observed with zero SiO2 interlayers. When the number of interlayers is incremented, a similar mode loss (close to −5 dB/cm) may be obtained for three or more interlayers. These results demonstrate that a WG fabricated with three SiO2 layers performs well for single mode guiding radiation at wavelengths from 1.5 µm to about 3 µm, which represents an ideal working range for our purposes.
The effective length (Leff) was simulated in order to work in the regime where linear loss still dominates over non-linear loss. In that regime, TPA and TPA-induced free-carrier absorption (FCA) were not produced. The second process is because free-carriers that are excited to the conduction band could further absorb photons and be excited to higher energy states. These carriers are distributed along the length of the device over the FCA effective length, LFCA,eff: In this context, we have simulated the cut-off wavelength for different interlayer numbers, as shown in Figure 5a. In this context, we have simulated the cut-off wavelength for different interlayer numbers, as shown in Figure 5a. Furthermore, in Figure 5b, the mode loss and the effective index at cut-off wavelengths are plotted. A dramatic increase in mode loss is observed with zero SiO2 interlayers. When the number of interlayers is incremented, a similar mode loss (close to −5 dB/cm) may be obtained for three or more interlayers. These results demonstrate that a WG fabricated with three SiO2 layers performs well for single mode guiding radiation at wavelengths from 1.5 µm to about 3 µm, which represents an ideal working range for our purposes.
The effective length (Leff) was simulated in order to work in the regime where linear loss still dominates over non-linear loss. In that regime, TPA and TPA-induced free-carrier absorption (FCA) were not produced. The second process is because free-carriers that are excited to the conduction band could further absorb photons and be excited to higher energy states. These carriers are distributed along the length of the device over the FCA effective length, LFCA,eff: Furthermore, in Figure 5b, the mode loss and the effective index at cut-off wavelengths are plotted. A dramatic increase in mode loss is observed with zero SiO 2 interlayers. When the number of interlayers is incremented, a similar mode loss (close to −5 dB/cm) may be obtained for three or more interlayers. These results demonstrate that a WG fabricated with three SiO 2 layers performs well for single mode guiding radiation at wavelengths from 1.5 µm to about 3 µm, which represents an ideal working range for our purposes.
The effective length (L eff ) was simulated in order to work in the regime where linear loss still dominates over non-linear loss. In that regime, TPA and TPA-induced free-carrier absorption (FCA) were not produced. The second process is because free-carriers that are excited to the conduction band could further absorb photons and be excited to higher energy states. These carriers are distributed along the length of the device over the FCA effective length, L FCA,eff : The factor of two in the denominator is due to the dependence of free-carrier density on the square of the optical intensity.
As shown in Figures 3 and 4, an increment of n eff occurs when we increment the number of SiO 2 interlayers. In order to understand the minimal number of SiO 2 interlayers to obtain a wide range of guided wavelengths, the propagation constant was calculated as a function of SiN and SiO 2 thickness for different numbers of layers ( Figure 6). The propagation constant (also known as transmission constant) is the measure of the change in amplitude and phase per unit distance of the field vector. Calculations performed show an increment in the propagation constant given by the layer thickness increment (Figure 6a,b), whilst a reduction in the propagation constant is produced by incrementing the number of layers constituting the WG (Figure 6c).

Appl. Sci. 2021, 11, x FOR PEER REVIEW 7 of 14
The factor of two in the denominator is due to the dependence of free-carrier density on the square of the optical intensity.
As shown in Figures 3 and 4, an increment of neff occurs when we increment the number of SiO2 interlayers. In order to understand the minimal number of SiO2 interlayers to obtain a wide range of guided wavelengths, the propagation constant was calculated as a function of SiN and SiO2 thickness for different numbers of layers ( Figure 6). The propagation constant (also known as transmission constant) is the measure of the change in amplitude and phase per unit distance of the field vector. Calculations performed show an increment in the propagation constant given by the layer thickness increment ( Figure  6a,b), whilst a reduction in the propagation constant is produced by incrementing the number of layers constituting the WG (Figure 6c).  Figure 7a shows the effective length (Leff) calculated for the cut-off refractive index for different numbers of interlayers. We can observe that Leff differs from the linear behavior, (i.e., the FCA losses are higher). In particular, in Figure 7b, Leff for three SiO2 interlayers is plotted. Leff of the device is a function of the physical length (L) and the linear loss (αlin). In that regime, the linear loss still dominates over the non-linear loss and, thus, we can calculate the optimal WG length. We can note that above 6 mm, Leff starts to draw away from the linear behavior. That value for Leff is similar to that shown by Motamedi et al. [30] (Leff = 6 mm) for Si WGs. In their study, Si WGs had a dimension of 106 nm × 475 nm with  Figure 7a shows the effective length (L eff ) calculated for the cut-off refractive index for different numbers of interlayers. We can observe that L eff differs from the linear behavior, (i.e., the FCA losses are higher). In particular, in Figure 7b, L eff for three SiO 2 interlayers is plotted. L eff of the device is a function of the physical length (L) and the linear loss (α lin ). In that regime, the linear loss still dominates over the non-linear loss and, thus, we can calculate the optimal WG length. We can note that above 6 mm, L eff starts to draw away from the linear behavior. That value for L eff is similar to that shown by Motamedi et al. [30] (L eff = 6 mm) for Si WGs. In their study, Si WGs had a dimension of 106 nm × 475 nm with a length equal to 14.9 mm and guided a fundamental TE mode at the 1500 nm wavelength. That result indicates that 6 mm is a good WG length for our purposes.  Once we established the number of layers, we calculated the range of WG widths that offer less loss, as is shown in Figure 8. In Figure 8b, losses are plotted as a function of wavelength and WG width. In particular, for WG widths below 6 µm, a reduction in the range of guided wavelengths is obtained. As is shown in Figures 5a and 8a, for WGs fabricated with three SiO2 interlayers, the cut-off wavelength is about 3.1 µm. This particular can be solved for WG widths wider than 6 µm (Figure 8b) to obtain a wider range of guided wavelengths (from 300 nm to 3.1 µm). In order to perform an experimental test of the waveguiding behavior by using laser sources emitting in the visible and telecom wavelengths, we requested of the manufacturer five WGs with widths ranging between 7 and 10 µm.
For the 9.4 µm-width WG, we report in Figure 9 the simulation (left) of fundamental TE/TM modes for 1550 nm and 2640 nm. Once we established the number of layers, we calculated the range of WG widths that offer less loss, as is shown in Figure 8. In Figure 8b, losses are plotted as a function of wavelength and WG width. In particular, for WG widths below 6 µm, a reduction in the range of guided wavelengths is obtained. As is shown in Figures 5a and 8a, for WGs fabricated with three SiO 2 interlayers, the cut-off wavelength is about 3.1 µm. This particular can be solved for WG widths wider than 6 µm (Figure 8b) to obtain a wider range of guided wavelengths (from 300 nm to 3.1 µm). a length equal to 14.9 mm and guided a fundamental TE mode at the 1500 nm wavelength. That result indicates that 6 mm is a good WG length for our purposes. Once we established the number of layers, we calculated the range of WG widths that offer less loss, as is shown in Figure 8. In Figure 8b, losses are plotted as a function of wavelength and WG width. In particular, for WG widths below 6 µm, a reduction in the range of guided wavelengths is obtained. As is shown in Figures 5a and 8a, for WGs fabricated with three SiO2 interlayers, the cut-off wavelength is about 3.1 µm. This particular can be solved for WG widths wider than 6 µm (Figure 8b) to obtain a wider range of guided wavelengths (from 300 nm to 3.1 µm). In order to perform an experimental test of the waveguiding behavior by using laser sources emitting in the visible and telecom wavelengths, we requested of the manufacturer five WGs with widths ranging between 7 and 10 µm.
For the 9.4 µm-width WG, we report in Figure 9 the simulation (left) of fundamental TE/TM modes for 1550 nm and 2640 nm. In order to perform an experimental test of the waveguiding behavior by using laser sources emitting in the visible and telecom wavelengths, we requested of the manufacturer five WGs with widths ranging between 7 and 10 µm.
For the 9.4 µm-width WG, we report in Figure 9 the simulation (left) of fundamental TE/TM modes for 1550 nm and 2640 nm.  Figure (b,d) show the modes at 1550 nm and 2640 nm, respectively, obtained with an InSb camera for near infrared (NIR)/mid infrared (MIR) light.
A thorough WG characterization was performed by using the four lasers previously described. The inset of Figure 2 shows the WG behavior in the visible range, at 635 nm. All the light observed in the picture is solely due to a scattering process along the WG. Figure 9 (right) shows the NIR and MIR beam profiles guided through the WG, as recorded by the InSb camera. Such images, in agreement with the simulations (Figure 9 left), are a demonstration of single-mode waveguiding in the IR.

Finite Element Modeling for Stress Characterization
In order to get information on the local second-order non-linearity, a study of the behavior and influence of stress is necessary. Again, a stress evaluation across the multilayer WG cross section was made by using a FEM approach. We started this analysis from the previous analysis of Section 3.1.1.
The interlayer surface was defined with a surface load obtained according to the mass over the layer. A finer, free-triangular grid was employed to mesh the computational model. Deformation, stress, and strain at static equilibrium, in a stationary regime, were calculated assuming the following equations [59]: The total stress (σ) is given by: where the additional stress (σad) is the sum of initial stress (σ0), external stress (σext), and q stress (σq): The total strain contribution (ε) is established by the elastic (εel) and inelastic (εinel) strains. Inelastic strain (εinel) is given by the sum of initial strain (ε0), thermal strain (εth), hygroscopic strain (εhs), plastic strain (εpl), and creep strain (εcp) [59]. A thorough WG characterization was performed by using the four lasers previously described. The inset of Figure 2 shows the WG behavior in the visible range, at 635 nm. All the light observed in the picture is solely due to a scattering process along the WG. Figure 9 (right) shows the NIR and MIR beam profiles guided through the WG, as recorded by the InSb camera. Such images, in agreement with the simulations (Figure 9 left), are a demonstration of single-mode waveguiding in the IR.

Finite Element Modeling for Stress Characterization
In order to get information on the local second-order non-linearity, a study of the behavior and influence of stress is necessary. Again, a stress evaluation across the multilayer WG cross section was made by using a FEM approach. We started this analysis from the previous analysis of Section 3.1.1.
The interlayer surface was defined with a surface load obtained according to the mass over the layer. A finer, free-triangular grid was employed to mesh the computational model. Deformation, stress, and strain at static equilibrium, in a stationary regime, were calculated assuming the following equations [59]: The total stress (σ) is given by: where the additional stress (σ ad ) is the sum of initial stress (σ 0 ), external stress (σ ext ), and q stress (σ q ): The total strain contribution (ε) is established by the elastic (ε el ) and inelastic (ε inel ) strains. Inelastic strain (ε inel ) is given by the sum of initial strain (ε 0 ), thermal strain (ε th ), hygroscopic strain (ε hs ), plastic strain (ε pl ), and creep strain (ε cp ) [59]. Figure 10 shows the color maps for the stress components σ xx and σ yy obtained from the equations. Several interesting features can be seen from the maps. In order to evaluate the stress-induced behavior, we simulated the stress profile along the x-and y-coordinates. In Figure 10c, it can be easily observed that the σ xx component is not uniform along the y-coordinate; a critical increment is detected for SIN y layers with relaxation behavior in SiO 2 layers introducing a stress gradient.
Appl. Sci. 2021, 11, x FOR PEER REVIEW 10 of 14 Figure 10 shows the color maps for the stress components σxx and σyy obtained from the equations. Several interesting features can be seen from the maps. In order to evaluate the stress-induced behavior, we simulated the stress profile along the x-and y-coordinates. In Figure 10c, it can be easily observed that the σxx component is not uniform along the y-coordinate; a critical increment is detected for SINy layers with relaxation behavior in SiO2 layers introducing a stress gradient.

Measurements
Using the setup described in Section 2.2, we realized a proof of principle test of WG loss in terms of transmitted power of some of the modeled structures on a large spectral range (from visible to MIR). In Figure 11, we report the average losses vs. the WG width of 6-mm-length SiN multilayer structures,for four different wavelengths. These values indicate that losses have a strong dependence on the wavelength of the propagating radiation and have a quasi-non-dependence on that range of WG widths, as expected. Moreover, from the measurements and from our FEM analysis, such WGs are demonstrated to work simultaneously in the NIR/MIR in single spatial mode.

Measurements
Using the setup described in Section 2.2, we realized a proof of principle test of WG loss in terms of transmitted power of some of the modeled structures on a large spectral range (from visible to MIR). In Figure 11, we report the average losses vs. the WG width of 6-mm-length SiN multilayer structures, for four different wavelengths. These values indicate that losses have a strong dependence on the wavelength of the propagating radiation and have a quasi-non-dependence on that range of WG widths, as expected. Moreover, from the measurements and from our FEM analysis, such WGs are demonstrated to work simultaneously in the NIR/MIR in single spatial mode.
Appl. Sci. 2021, 11, x FOR PEER REVIEW 10 of 14 Figure 10 shows the color maps for the stress components σxx and σyy obtained from the equations. Several interesting features can be seen from the maps. In order to evaluate the stress-induced behavior, we simulated the stress profile along the x-and y-coordinates. In Figure 10c, it can be easily observed that the σxx component is not uniform along the y-coordinate; a critical increment is detected for SINy layers with relaxation behavior in SiO2 layers introducing a stress gradient.

Measurements
Using the setup described in Section 2.2, we realized a proof of principle test of WG loss in terms of transmitted power of some of the modeled structures on a large spectral range (from visible to MIR). In Figure 11, we report the average losses vs. the WG width of 6-mm-length SiN multilayer structures,for four different wavelengths. These values indicate that losses have a strong dependence on the wavelength of the propagating radiation and have a quasi-non-dependence on that range of WG widths, as expected. Moreover, from the measurements and from our FEM analysis, such WGs are demonstrated to work simultaneously in the NIR/MIR in single spatial mode. Figure 11. Coupling losses measured at 532 nm, 635 nm, 1550 nm, and 2640 nm for each WG width. Figure 11. Coupling losses measured at 532 nm, 635 nm, 1550 nm, and 2640 nm for each WG width.

Discussion
In the IR range, losses are below 6 dB (measured at 1550 and 2640 nm), while, in the visible range, we observe an increment of about 20 dB (measured at 532 and 635 nm). Insertion losses are the sum of two contributions: propagation losses, probably due to different scattering processes along the WGs; coupling losses, due to spatial mode mismatch; and Fresnel reflections at input/output facets. Coupling losses cannot be easily estimated, due to the complexity of the WGs' geometry, preventing reliable numerical calculations. Nevertheless, we can set an upper limit for the coupling losses, assuming them to be equal to the coupling losses at 2640 nm, since at this longer wavelength, propagation losses are expected to be the lowest ones. As we indicate from the results obtained in Section 3.2, the lengths of the guides used in the test were below the effective length (6 mm), and, therefore, we assume that the propagation losses were negligible, and the losses measured are dueto coupling loss only. Losses measured at 1550 nm were close to those reported in other studies, as such 5.6 dB for Si WGs at 1500 nm with a 106 nm × 475 nm core [27] or 8 dB for SiN WGs at 1540 nm with a cross-section of 1.11 × 1.11 µm 2 [28].
Using FEM, we simulated different lengths for every WG and, consequently, we evaluated the propagation constant. In these calculations, the number of layers, the width, or the length could be fixed. In particular, considering 6-mm guides with three SiO 2 interlayers of different widths, a good agreement in terms of order of magnitude was obtained between simulation ( Figure 8) and experimental evaluation of propagation losses ( Figure 11). The measured values show a loss reduction at longer wavelengths. Furthermore, a slight increase in transmission for wider WGs was observed. Such increments are associated with better light coupling, due to the wider lateral dimensions.
Furthermore, as can be seen in Figures 7 and 11, this particular geometry, optimized for single-mode waveguiding between 1 and 3 µm, shows low losses and single-mode behavior also at shorter wavelengths, around 1.5 µm. This is due to a trade-off between the width of the guide, the number and thickness of the multilayers, and the scattering losses. The FEM simulations proved to be a very powerful tool to design this kind of photonic structure. Indeed, a careful modelling of multi-wavelength guides will enable modal phase-matching towards more efficient non-linear processes or for integrated photonic circuits useful in quantum optics, communication, spectroscopy, and sensing experiments operating simultaneously among NIR and MIR regions.

Conclusions
In this study, detailed theoretical modeling and experimental measurements were shown for novel silicon-based WG geometry. SiN multilayer structures were fabricated based on finite element modeling. At first, we fixed a set of parameters (wavelength and polarization) and performed simulations by varying interlayer height and WG width. The simulations show efficient single-mode waveguiding at wavelengths in the IR range (1.5-3 µm) for a WG made of 125-nm-thick SiN layers and 90-nm-thick SiO 2 layers. Furthermore, a stress simulation was performed, showing a stress gradient in the center of every WG. We characterized the waveguiding behavior in a wide range of wavelengths (532, 635, 1550, and 2640 nm). The WGs showed good waveguiding behavior at NIR-MIR wavelengths (both 1550 and 2640 nm, with loss values close to 4 dB), while losses increased at wavelengths in the visible range (about 20 dB). The losses measured for multilayer WGs were similar to those reported by other studies for Si or Si 3 N 4 WGs [27,28]. Instead, from calculations, we did not observe any variation in loss values for different WG widths, because of minimal variations in the stress profile.
In addition, ZBLAN and chalcogenide glass are emerging as materials for MIR, particularly for fiber components. However, silicon technologies enable the realization of planar photonic circuits and devices (i.e., filters or resonators). The coupling, not simple in such a spectral region, could be improved by using precisely designed chalcogenide optical components.
In summary, we can conclude that Si x N y /SiO 2 multilayer WGs are good candidates for simultaneous MIR/NIR low-cost integrated optical components, providing broadband operation for spectroscopy and quantum technologies.
Author Contributions: For research articles I.L.G. performed the measurements, performed the numerical simulations, analyzed the data, and wrote the manuscript; M.S.d.C. contributed to the conception of the experiment, to the development of the experimental setup, and to simulation analysis; I.G. performed experimental measurements; P.C.P. and D.M. contributed to the conception of the experiment, to the simulations, and to the data interpretation; and P.D.N. contributed to the conception of the experiment and coordinated the experiment. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement:
The data presented in this study are available on request from the authors.