Issue 
J. Eur. Opt. SocietyRapid Publ.
Volume 20, Number 1, 2024
Current challenges and solutions in guided wave optics and integrated photonics



Article Number  9  
Number of page(s)  9  
DOI  https://doi.org/10.1051/jeos/2024007  
Published online  01 April 2024 
Research Article
Multiple Bloch surface wave excitation with gratings
^{1}
Center for Photonics Sciences, Department of Physics and Mathematics, University of Eastern Finland, PO Box 111, 80101 Joensuu, Finland
^{2}
Dipartimento di Scienza Applicata e Tecnologia, Politecnico di Torino, 10129 Torino, Italy
^{3}
Dispelix Oy, Yliopistokatu 7, FI80130 Joensuu, Finland
^{*} Corresponding author: lewisa@uef.fi
Received:
14
December
2023
Accepted:
21
February
2024
We study the coupling of a finite number of Bloch Surface Waves (BSWs) propagating in different directions at the surface of a dielectric multilayer. These surface waves arise from a set of diffraction orders associated to a grating on the bottom surface of the substrate that is illuminated by a normally incident beam. Simultaneous excitation of multiple BSWs is possible with a set of diffraction orders having the same radial spatial frequency. Using rigorous electromagnetic theory, we design gratings for simultaneous excitation of two, four and six BSWs propagating in directions separated by π, π/2 and π/3 azimuthal intervals, respectively.
Key words: Evanescent waves / Surface electromagnetic waves / Bloch surface waves / Multiple bloch surface waves / MBSWs
© The Author(s), published by EDP Sciences, 2024
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1 Introduction
Surface Electromagnetic Waves (SEW) represent an interesting option for controlling optical signals on miniaturized chips for integrated optics and sensing applications. Surface Plasmon Polaritons (SPP) are probably the most widely known SEWs, but they exhibit inherent issues related to the ohmic losses introduced by the metallic materials involved. As an alternative, SEWs sustained by dielectric multilayers (ML) have attracted a growing interest in the past decade. This kind of SEW [1] is also referred to as Bloch Surface Waves (BSWs) [2] to highlight the role of the underlying periodic multilayer structure required for their existence. BSWs offer several advantages as compared to SPPs, such as a wide spectral tunability and low losses thanks to the large choice of transparent dielectric materials available for multilayer manufacturing. In addition, BSWs can be either TE or TMpolarized [3, 4], depending on the multilayer design. Their excitation by pulsed fields has also been recently studied numerically [5].
Being surface waves, BSWs are evanescent in the medium above the multilayer surface. The coupling with freespace radiation in a BSWbased device is therefore critical as it must provide momentum matching beyond the lightline. In most of the applications proposed so far, BSW coupling is performed by means of bulky prisms, either in Kretschmann or Otto configuration [6]. However, more sophisticated approaches have been recently implemented, involving, for example, the use of individual scatterers [7–9] or miniature prisms [10] placed onto the multilayer surface. Another promising option is represented by integrating diffraction gratings within the BSWsupporting structure [11]. This has been done mainly in two different ways: with the grating being fabricated on top of the multilayer [12–14] or buried beneath the multilayer [15]. In the first case, the multilayer is substantially planar, with the exception of the top layer, where the grating unavoidably perturbs the dispersion of the BSW mode (dielectric loading/unloading effect). In the second case, the grating is fabricated on the substrate surface prior to the multilayer deposition, which occurs on the same side. The multilayer itself results to not be perfectly planar because it (partially) conforms to the underlying corrugation. In both configurations, the BSW dispersion is altered by the presence of the grating, which may lead to some difficulties regarding the precise control of optical functions of complex, possibly resonating, BSWbased architectures.
We propose an alternative approach on diffractive coupling for BSWs, with gratings fabricated on the bottom surface of a transparent substrate having the multilayer deposited on the top surface. In particular, we explore the possibility of using twodimensional gratings to simultaneously couple BSWs propagating in more than two directions by exploiting the momentum distribution of several diffraction orders. Once the mode dispersion of the multilayer is known, our approach facilitates BSW coupling in a controllable way, as far as wavelengths/numbers and propagation directions are concerned. The directional coupling of BSWs has been already tackled in a few previous articles [16–18], although never considered for multiple directions at once. When the optical path through the substrate is also taken into account, our approach allows a predictable control onto the coupling locations of BSWs launched in different directions.
The present paper is composed as follows. We begin, in Section 2, by introducing the gratingbased BSW excitation principle and the assumed geometrical configuration. The theoretical framework for grating design, for which we use a rigorous technique known as the Fourier Modal Method (FMM) [19], is described in Section 3. The design process is analogous with the synthesis of gratingbased multiple freespace beam splitters [20], but here we need to account for the BSW excitation conditions and the polarization state of the input wave. In Section 4, we first consider BSW stack design, providing a “benchmark” stack employed in the rest of the work, and then cover the design of linear gratings for simultaneous excitation of two counterpropagating BSWs. Such designs are extended in Section 5 to twodimensional periodic gratings for excitation of either four or six BSWs propagating at 90° or 60° intervals along the stack, respectively. After a discussion presented in Section 6, conclusions are drawn in Section 7.
2 Excitation principle and geometry
Figure 1 illustrates the geometry for the simplest case of excitation of two counterpropagating BSWs. A flat fused silica substrate with refractive index n_{sub} = 1.462, such as a 0.5 mm or 3 mmthick SiO_{2} plate, is illuminated by a normally incident monochromatic beam (wavelength λ_{0}) from the medium underneath (air). A linear grating, with period d of the order of λ_{0}, provided on the airsubstrate surface, splits the beam into three transmitted orders propagating within the substrate: the zeroth order m = 0 and the first diffracted orders m = ±1. The orders m = ±1 propagate in directions θ_{+1} and θ_{−1} given by sinθ_{±1} = ±λ_{0}/n_{sub}d towards the multilayer stack on the top surface of the substrate. If θ_{±1} matches the Kretschmannincidence BSW excitation angle θ_{BSW} for the given wavelength and polarization state (TE or TM), two counterpropagating BSWs are generated simultaneously. The excitation is efficient as long as the angular spectrum of each diffracted order, which defines the beam divergence, falls essentially within the (stackdependent) BSW momentum bandwidth. The polarization state of illumination affects the coupling significantly; we will consider only BSW excitation in TE polarization, which generally requires a smaller number of stack layers than TMpolarized BSW excitation.
Figure 1 Principle of MBSW generation: the twobeam case. A binary linear surfacerelief grating defined by period d, ridge width c, and ridge height h on the bottom of a substrate of thickness H splits the input beam into two diffracted orders m = −1 and m = +1, which excite BSWs on the top surface of the substrate by interaction with the multilayer stack (ML). We assumed n_{sub} = 1.462 and n_{sup} is air. 
The parameters of the system are chosen such that the two BSWs shown in Figure 1 are spatially separated under finitebeam illumination. This feature can be useful in BSWbased platforms such as interferometers [21] and integrated components [22]. Firstorder diffracted beams are partially reflected at the top interface, thus propagating back into the substrate. The reflected beams continue to propagate according to multiplereflection paths inside the substrate unless they are extracted by means of diffusers or gratings. At each reflection with the ML interface, coupling to BSW occurs. Stated differently, BSWs are launched at different locations on the ML surface each time the beam is incident on the bottom interface of the dielectric stack, thus leading to the appearance of BSW interference effects unless the substrate thickness H is sufficiently large to minimize spatial overlaps.
3 Theoretical framework
Let us consider a rectangularly periodic grating of period d_{x} × d_{y} in the cartesian xy coordinate system and assume a planewave illumination (at frequency ω) normally incident onto the substrate from air. In view of the grating equations, the wave vectors of the propagating diffraction orders (m, n) in the substrate are$${\mathit{k}}_{\mathrm{mn}}={k}_{\mathrm{xm}}\widehat{\mathit{x}}+{k}_{\mathrm{yn}}\widehat{\mathit{y}}+{k}_{\mathrm{zmn}}\widehat{\mathit{z}}$$(1)where$${k}_{\mathrm{xm}}=m{K}_{x}=2\pi m/{d}_{x}$$(2a) $${k}_{\mathrm{yn}}=n{K}_{y}=2\pi n/{d}_{y}$$(2b) $${k}_{\mathrm{zmn}}=\sqrt{{k}_{0}^{2}{n}_{\mathrm{sub}}^{2}{k}_{\mathrm{xm}}^{2}{k}_{\mathrm{yn}}^{2}}$$(2c)and k_{0} = ω/c_{0} = 2π/λ_{0} is the wave number in vacuum. After defining the radial spatial frequency of the generic order (m, n) as$${k}_{\rho \mathrm{mn}}=\sqrt{{k}_{\mathrm{xm}}^{2}+{k}_{\mathrm{yn}}^{2}}=\sqrt{(m{K}_{x}{)}^{2}+(n{K}_{y}{)}^{2}}$$(3)the condition k_{ρmn} < k_{0}n_{sub} identifies those diffraction orders propagating within the substrate, the others being evanescent. If we denote the refractive index of the superstrate by n_{sup} and assume n_{sup} < n_{sub}, order (m, n) is evanescent in the superstrate when k_{ρmn} > k_{0}n_{sup}. Considering BSW excitation, we are therefore interested in orders with radial spatial frequencies in the range k_{0}n_{sup} < k_{ρmn} < k_{0}n_{sub}. We are primarily interested in the nearest neighbors of the zeroth transmitted order, while higher orders are made evanescent by appropriate choices of d_{x} and d_{y}. In the illustrative example presented in Figure 2a, orders (m, n) = (−1, 0) and (m, n) = (+1, 0) fall on the yellow line of radius k_{ρBSW}, which defines the BSW excitation condition dictated by the ML design.
Figure 2 (a) Diffraction orders of a rectangular lattice in spatialfrequency representation at normal incidence. Diffraction orders are represented by dots at positions k_{xm} = mK_{x}, k_{y} = nK_{y}. Blue and red circles represent the cutoff radial spatial frequencies k_{ρ} = k_{0}n_{sub} and k_{ρ} = k_{0}n_{sup}, respectively, between which BSW excitation is possible. The yellow circle indicates the radial spatial frequency of BSW for a given ML. (b) Definition of the propagation angles (θ_{mn}, ϕ_{mn}) of a single transmitted diffracted order (m, n) in its exit plane (the grey rectangle) and the π − σ basis of the diffracted electric field. 
Following reference [23], we define the “exit plane” of diffraction order (m, n) as the plane containing the wave vector k _{mn} and the unit vector $\widehat{\mathit{z}}$. Further, propagation angles θ_{mn} and ϕ_{mn} of the transmitted orders, are defined as$${k}_{\mathrm{xm}}={k}_{0}{n}_{\mathrm{sub}}\mathrm{sin}{\theta}_{\mathrm{mn}}\mathrm{cos}{\varphi}_{\mathrm{mn}}$$(4a) $${k}_{\mathrm{yn}}={k}_{0}{n}_{\mathrm{sub}}\mathrm{sin}{\theta}_{\mathrm{mn}}\mathrm{sin}{\varphi}_{\mathrm{mn}}$$(4b) $${k}_{\mathrm{zmn}}={k}_{0}{n}_{\mathrm{sub}}\mathrm{cos}{\theta}_{\mathrm{mn}}$$(4c)as illustrated in Figure 2b. Here ϕ_{mn} is the azimuthal angle in the range [0, 2π), measured counterclockwise from the k_{x} axis, and θ_{mn} in the range [0, π/2) is the propagation angle measured from the k_{z} axis. It will prove convenient to use the socalled π–σ basis (or local TM/TE basis) to define the polarization states of the transmitted orders. As described in reference [23], this basis allows us to treat incident fields with any polarization state, including partial polarization. Here, however, we are mainly interested in either fully polarized or unpolarized illumination.
If the incident plane wave is fully polarized, we can use any suitable rigorous grating analysis method (in our case FMM) to determine the transverse Cartesian components e_{xmn} and e_{ymn} of the polarization vector for any transmitted order, as discussed shortly below. The longitudinal component of e _{mn} is fixed by Maxwell’s divergence equation, which gives k _{mn} ∙ e _{mn} = 0 and$${e}_{\mathrm{zmn}}=\frac{1}{{k}_{\mathrm{zmn}}}\left({k}_{\mathrm{xm}}{e}_{\mathrm{xmn}}+{k}_{\mathrm{yn}}{e}_{\mathrm{ymn}}\right).$$(5)
In the π–σ basis the polarization state of any order is described by a twodimensional vector e _{πσmn} = [e_{πmn}, e_{σmn}]^{T}, where the π and σ components are explicitly given by$${e}_{\pi \mathrm{mn}}={e}_{\mathrm{xmn}}\mathrm{cos}{\theta}_{\mathrm{mn}}\mathrm{cos}{\varphi}_{\mathrm{mn}}+{e}_{\mathrm{ymn}}\mathrm{cos}{\theta}_{\mathrm{mn}}\mathrm{sin}{\varphi}_{\mathrm{mn}}{e}_{\mathrm{zmn}}\mathrm{sin}{\theta}_{\mathrm{mn}}$$(6a) $${e}_{\sigma \mathrm{mn}}={e}_{\mathrm{xmn}}\mathrm{sin}{\varphi}_{\mathrm{mn}}+{e}_{\mathrm{ymn}}\mathrm{cos}{\varphi}_{\mathrm{mn}}.$$(6b)
As shown in Figure 2b, the component e_{πmn} lies in the exit plane, whereas e_{σmn} is perpendicular to it. Hence, they represent the TM and TE components of the electric field in the exit plane, respectively.
In diffraction by twodimensionally periodic gratings, the polarization states of the transmitted (and reflected) diffracted orders generally depend on the state of input polarization. We represent the polarization vector of a (generally, elliptically polarized) unitamplitude input plane wave as$$\mathit{e}={e}_{x}\widehat{\mathit{x}}+{e}_{y}\widehat{\mathit{y}}=\widehat{\mathit{x}}\mathrm{cos}\alpha +\widehat{\mathit{y}}\mathrm{sin}\alpha \mathrm{exp}\left(\mathrm{i\delta}\right)$$(7)normalized such that e = 1. The effect of the grating on transmitted radiation can be analyzed by calculating (by FMM) the transmission coefficients$${T}_{\mathrm{xmn}}^{\left(x\right)},{T}_{\mathrm{ymn}}^{\left(x\right)},{T}_{\mathrm{xmn}}^{\left(y\right)},{T}_{\mathrm{ymn}}^{\left(y\right)}$$(8)
for all diffraction orders, where the superscripts (x) and (y) refer to illumination by a purely xpolarized ($\mathit{e}=\widehat{\mathit{x}}$) or ypolarized ($\mathit{e}=\widehat{\mathit{y}}$) incident wave. The coefficients in equation (8) are precisely the complex vector amplitudes that appear in the Rayleigh planewave expansion of the field at the output plane of the grating; see, e.g., equation (5) in reference [20]. For an arbitrarily (fully) polarized incident wave the transverse electricfield components of the transmitted orders are [23]$${e}_{\mathrm{xmn}}={T}_{\mathrm{xmn}}^{\left(x\right)}{e}_{x}+{T}_{\mathrm{xmn}}^{\left(y\right)}{e}_{y}$$(9a) $${e}_{\mathrm{ymn}}={T}_{\mathrm{ymn}}^{\left(x\right)}{e}_{x}+{T}_{\mathrm{ymn}}^{\left(y\right)}{e}_{y}.$$(9b)
The longitudinal components e_{zmn} are obtained from equation (5), and the π–σ representation of each order is given by equations (6). Since the input polarization state affects both the π and σ components, it can be used as a design degree of freedom in multipleBSW excitation, in addition to the geometrical grating parameters.
It is customary to describe the state of polarization of a fully polarized field by a 2 × 2 polarization matrix J = e ^{*} e ^{T} ([24], Sec. 6.3.2). Explicitly, for the incident field,$$\mathbf{J}=\left[\begin{array}{ll}{J}_{\mathrm{xx}}& {J}_{\mathrm{xy}}\\ {J}_{\mathrm{yx}}& {J}_{\mathrm{yy}}\end{array}\right]=\left[\begin{array}{ll}{\left{e}_{x}\right}^{2}& {e}_{x}^{\mathrm{*}}{e}_{y}\\ {e}_{y}^{\mathrm{*}}{e}_{x}& {\left{e}_{y}\right}^{2}\end{array}\right]$$(10)where the asterisk denotes complex conjugation. Correspondingly, the polarization state of any transmitted order in the π–σ basis is described by ${\mathbf{J}}_{\pi \sigma \mathrm{mn}}={\mathbf{e}}_{\pi \mathrm{mn}}^{*}{\mathbf{e}}_{\sigma \mathrm{mn}}^{\mathrm{T}}$ [23], explicitly$${\mathbf{J}}_{\pi \sigma \mathrm{mn}}=\left[\begin{array}{ll}{J}_{\pi \pi \mathrm{mn}}& {J}_{\pi \sigma \mathrm{mn}}\\ {J}_{\sigma \pi \mathrm{mn}}& {J}_{\sigma \sigma \mathrm{mn}}\end{array}\right]=\left[\begin{array}{ll}{\left{e}_{\pi \mathrm{mn}}\right}^{2}& {e}_{\pi \mathrm{mn}}^{*}{e}_{\sigma \mathrm{mn}}\\ {e}_{\sigma \mathrm{mn}}^{*}{e}_{\pi \mathrm{mn}}& {\left{e}_{\sigma \mathrm{mn}}\right}^{2}\end{array}\right].$$(11)
The polarization states of the diffracted orders can also be characterized by the Stokes parameters [23]$${S}_{0\mathrm{mn}}={J}_{\pi \pi \mathrm{mn}}+{J}_{\sigma \sigma \mathrm{mn}}$$(12a) $${S}_{1\mathrm{mn}}={J}_{\pi \pi \mathrm{mn}}{J}_{\sigma \sigma \mathrm{mn}}$$(12b) $${S}_{2\mathrm{mn}}=2\mathfrak{R}\left({J}_{\pi \sigma \mathrm{mn}}\right)$$(12c) $${S}_{3\mathrm{mn}}=2\mathfrak{I}\left({J}_{\pi \sigma \mathrm{mn}}\right)$$(12d)where $\mathfrak{R}$ and $\mathfrak{I}$ denote the real and imaginary parts. The normalized forms of the Stokes parameters are defined as s_{jmn} = S_{jmn}/S_{0mn} (j = 1, 2, 3), and the degree of polarization associated with order (m, n) is given by$${P}_{\mathrm{mn}}=\sqrt{{s}_{1\mathrm{mn}}^{2}+{s}_{2\mathrm{mn}}^{2}+{s}_{3\mathrm{mn}}^{2}}.$$(13)
For a fully polarized incident wave, P_{mn} = 1 for all orders, even though the values of the individual Stokes parameters generally depend on order indices.
In addition to fully polarized illumination, we consider the opposite extreme case of unpolarized illumination. The matrix J for partially polarized light is defined as J = 〈e ^{*} e ^{T}〉, where the brackets denote ensemble averaging over all polarization realizations. For unpolarized illumination it has a diagonal form ([24], Sec. 6.3.3)$$\mathbf{J}=\frac{1}{2}\left[\begin{array}{ll}1& 0\\ 0& 1\end{array}\right]$$(14)and the degree of input polarization is P = 0. The polarization matrix associated with order (m, n) can be represented as an average$${\mathbf{J}}_{\pi \sigma \mathrm{mn}}=\frac{1}{2}\left[{\mathbf{J}}_{\pi \sigma \mathrm{mn}}^{\left(x\right)}+{\mathbf{J}}_{\pi \sigma \mathrm{mn}}^{\left(y\right)}\right]$$(15)which remains diagonal because e_{x} and e_{y} are uncorrelated. However, since the grating treats these components differently, in general J_{ππmn} ≠ J_{σσmn}, implying that the individual orders become partially polarized with P_{mn} > 0.
In standard beam splitting problems in resonancedomain diffractive optics [20] one is interested in the distribution of the diffraction efficiencies of the propagating orders. Since we have normalized the intensity of the incident field such that S_{0} = 1, the diffraction efficiencies are defined as [23]$${\eta}_{\mathrm{mn}}={n}_{\mathrm{sub}}\mathrm{cos}{\theta}_{\mathrm{mn}}{S}_{0\mathrm{mn}}.$$(16)
In BSW excitation problem, the design goal is to maximize and equalize the coupling of the incident field to a set of BSW modes with the angle θ_{mn} equal to θ_{BSW}. If θ_{BSW} is the excitation angle for TE polarization, the component e_{σmn} excites a BSW while e_{πmn} is nonresonant, and vice versa.
We choose the geometry such that several diffraction orders have the same radial spatial frequency k_{ρmn} = k_{0}n_{sub}sinθ_{BSW}, and therefore lie on the yellow circle depicted in Figure 2a. The relative amplitudes of the excited BSWs are determined by the σpolarized components J_{σσmn} in TE polarization and J_{ππmn} in TM polarization. The fraction$${\kappa}_{\mathrm{mn}}={J}_{\sigma \sigma \mathrm{mn}}/{J}_{\pi \pi \mathrm{mn}}$$(17)provides the ratio of the coupled and uncoupled parts of the incident wave in BSW excitation.
4 Planewave design with linear gratings
As evident from the preceding discussion, the π–σ representation returns the BSW excitation problem to the basic TE or TM polarized problem. In addition, since we assume a substrate thickness H ≫ λ_{0}, the evanescent parts of the diffracted fields above the grating and the BSW field below the stack are spatially well separated. Hence, we may treat the BSW stack design and the grating design as two separate problems. In order to obtain an illustrative stack design useful for our purposes, we fix λ_{0} = 514 nm and provide a stack geometry sustaining BSWs at angles between the blue and red lines in Figure 2a. The resulting stack can then be used to design gratings for excitation of BSWs that lie on the yellow circle in Figure 2a.
4.1 Multilayer stack design
Figure 3a shows the assumed stack structure, which consists of N high/low (H/L) refractive index bilayers and a terminating top (T) layer with refractive indices n_{H}, n_{L}, n_{T} and thicknesses h_{H}, h_{L}, h_{T}, respectively. To reduce the number of variable parameters, we consider TE polarization, fix the number of bilayers to N = 6, use refractive indices n_{H} = 2.520 (TiO_{2}), n_{L} = 1.476 (SiO_{2}), n_{T} = n_{H} = 2.520. The thicknesses h_{H}, h_{L}, h_{T} are used to design the stack such that the BSW resonance occurs at an angle θ_{mn} in the exit plane.
Figure 3 (a) Definition of the multilayer structure and notation. (b) Stack parameters as a function of the ratio k_{x}/k_{0} = n_{eff} for TEmode BSW excitation with N = 6 bilayers: h_{H}/λ_{0} (blue), h_{L}/λ_{0} (red), and h_{T}/λ_{0} (black). The dots mark the position k_{x}/k_{0} = 1.1209 for BSW excitation at 50° angle of incidence. 
Figure 3b shows the design results. The horizontal axis is k_{x}/k_{0} = n_{sub} sinθ_{mn} = n_{eff}, where n_{eff} can be interpreted as the effective index of the stack. The plotting range starts from the critical angle of BSW generation and extends to k_{x}/k_{0}n_{sub}, i.e., it spans the region between the blue and red circles in Figure 2a. As n_{eff} increases, the BSW becomes increasingly buried within the multilayer and acts less like a surface mode. At the same time, all layer thicknesses show a monotonically increasing trend.
4.2 Twoway splitting
As illustrated in Figure 1, the coupling of two counterpropagating BSWs is possible with linear gratings (d_{x} = d, d_{y} = ∞). The exit plane of both orders, (m, n) = (−1, 0) and (+1, 0), is the xz plane and the π–σ representation reduces to the standard TM/TE decomposition. Since, by symmetry, η_{−1,0} = η_{+1,0} for binary profiles defined in the inset of Figure 1, we need to maximize η_{+1,0}. This also leads to the optimum value of J_{σσmn}, while J_{ππmn} = 0. Now we only need to find the values of the fill factor f = c/λ_{0} and grating height h/λ_{0} that maximize η_{+1,0} (= η_{−1,0}) to also maximize J_{σσ}.
The gratingdesign results are summarized in Figure 4. The optimum fill factor remains fairly constant over the entire angular range considered here, whereas the optimum grating height decreases with increasing angle. The efficiencies of all propagating orders are plotted in Figure 4b. At around k_{x}/k_{0} = 1.1209 (corresponding to an excitation angle 50°) we get η_{±1,0} ≈ 0.4973. Some light is “lost” in zeroth reflected and transmitted orders when we move close to the cutoff at k_{x}/k_{0} = 1 or towards larger values of k_{x}/k_{0}, but the designs remain acceptable over a relatively wide range of excitation angles.
Figure 4 Design of twoway beam splitters. (a) Optimum values of the fill factor f (red) and the relief depth h/λ_{0} (blue) as a function of the exit angle of the first diffracted order in TE polarization. (b) The corresponding firstorder diffraction efficiency η_{+1,0} = η_{−1,0} (black), the efficiency η_{0,0} of the zeroth transmitted order (red), and that of the zeroth reflected order (blue), which is 1 − 2η_{+1,0} − η_{0,0} due to energy conservation. The inset shows the grating structure and direction of illumination. 
The results in Figures 3b and 4 allow us to design twoway beam splitters for any BSW resonance angle of interest. The stack design for the desired angle is obtained from Figure 3b and the corresponding grating design from Figure 4a. The performance of the design can be evaluated from Figure 4b. To limit the number of variables further, we set θ_{BSW} = 50°, corresponding to k_{x}/k_{0} = 1.1209. The stack design is marked by the dots in Figure 3b, the optimum parameters for TE excitation with N = 6 bilayers being h_{H} = 60 nm, h_{L} = 85 nm, and h_{T} = 20 nm. Correspondingly, the vertical lines in Figure 4a give a grating design f = 0.2536, h/λ_{0} = 0.2752, with η_{±1,0} = 0.4973.
Considering the optimized case represented by the dots on Figure 3 and the black dashed line on Figure 4, simulation of the reflected and transmitted coefficients has been performed for the full structure. It implies a grating of period d ≃ 459 nm and fill factor is f = 0.2536 and h/λ_{0} = 0.2752 on the lower side of a fused silica wafer on top of which the multilayer is deposited. The multilayer design leads to a Bloch surface wave excited when the first diffracted order emerge from the grating at an angle of 50° (k_{x}/k_{0} = 1.1209) at a wavelength of 514 nm. This is observed in Figure 5a, where a strong dip in reflection arises at this value of k_{x}/k_{0}. In Figure 5b the response in wavelength is presented.
Figure 5 Response of the full structure (grating, substrate, multilayer and superstrate). (a) and (b) Reflected (black curves) and transmitted (red curves) first diffracted orders as a function of the normalized wavevector (a) and wavelength (b). 
5 Planewave design with biperiodic gratings
We proceed to design of twodimensionally periodic gratings that allow simultaneous excitation of more than two BSWs. Two lattice geometries are considered: square lattices for fourway excitation and hexagonal lattices for sixway excitation.
5.1 Fourway splitting
Let us first consider biperiodic gratings with primitive directlattice vectors ${\mathit{a}}_{1}=d\widehat{\mathit{x}}$, ${\mathit{a}}_{2}=d\widehat{\mathit{y}}$. The (Wigner–Seitz) primitive cell is squareshaped, covering the area −d/2 < x <d/2, −d/2 < y < d/2. The spatial frequencies of the diffraction orders are then k_{xm} = mK, k_{yn} = nK, the coordination number is 4, and the nearest neighbors of the zeroth order (0, 0), namely (m, n) = (+1, 0), (0, +1), (−1, 0), (0, −1), propagate in directions ϕ_{+1,0} = 0, ϕ_{0,+1} = π/2, ϕ_{−1,0} = π, ϕ_{0,−1} = 3π/2, respectively. By an appropriate choice of d, all of these four orders can be placed simultaneously on the yellow ring in Figure 2a, thus enabling fourway BSW excitation.
In the design, we found it sufficient to consider binary (zinvariant) relativepermittivity profiles of the particular form$${\u03f5}_{\mathrm{r}}\left(x,y,z\right)=\{\begin{array}{cc}{n}_{1}^{2}& \mathrm{when}{x}^{2}+{y}^{2}{r}^{2}\\ {n}_{2}^{2}& \mathrm{otherwise}\end{array}$$(18)in 0 < z < h within the primitive cell. The circular feature defined by the radius r can be either a pillar (n_{1} = n_{sub}, n_{2} = 1 or a hole (n_{1} = 1, n_{2} = n_{sub}) etched in the substrate. This type of pillar/hole structures can be patterned at a nanometerscale addressing resolution using electron beam lithography system available to us, and require only a single etching step.
The radius r and the relief depth h can be used as the structural design parameters. Some symmetry rules exist, which are helpful in the design. Since the unit cell and the structure are centered at the origin, the transmission coefficients in equation (8) satisfy the inversion symmetry rules$${T}_{x,m,n}={T}_{\mathrm{xmn}},{T}_{y,m,n}={T}_{\mathrm{ymn}}$$(19)for both (x) and (y) input polarizations. These rules hold regardless of the input polarization state, which however has an effect on the actual values of T_{xmn} and T_{ymn}. They reduce the number of orders that we need to (or can) control from four to two: we see from equation (19) that η_{−m,−n} = η_{mn}. Similar symmetry rules hold also for J_{σσmn} and J_{ππmn}.
We begin the design of fourway couplers by optimizing the structural parameters r and h to minimize the sum of the efficiencies of the reflected and transmitted zeroth orders. This maximizes the combined efficiency of the four nearestneighbor diffraction orders, and leaves the polarization state of the incident field free for design. Choosing θ_{BSW} = 50° (d ≈ 0.892λ_{0}), for either 45° or circularly polarized illumination and considering pillars, we get a design r ≈ 0.201λ_{0}, h ≈ 0.53176λ_{0}, which gives reflected and transmitted zeroorder efficiencies of ~3.5% and ~5.2%, respectively, leaving the rest of the incident energy to be distributed among the nearestneighbor orders.
Our remaining target is to equalize (and maximize) the coupling strengths J_{σσmn} of the four signal orders by designing the input polarization state defined in equation (7). The symmetry in the 4way splitting implies that there is no structurally induced polarization conversion: for (x)polarized input we get J_{σσmn} = 0 for orders (m, n) = (±1, 0), while (y)polarized input gives J_{σσmn} = 0 for orders (m, n) = (0, ±1). Considering linearly polarized light, the values of J_{σσmn} (and J_{ππmn}) vary rapidly with the angle α. Choosing α ≈ π/4 gives values J_{σσmn} ≈ 0.063 and κ_{mn} ≈ 0.359 for all four orders. The same result is obtained also for circularly polarized illumination with α ≈ π/4, δ = ±π/2. Both the optimized diffracted efficiencies and the maximized coupling strengths occur at the same illumination polarization.
Considering unpolarized illumination, the matrix J_{πσmn} becomes diagonal and the degree of polarization takes the form$${P}_{\mathrm{mn}}=\left{s}_{1\mathrm{mn}}\right=\frac{{J}_{\pi \pi \mathrm{mn}}{J}_{\sigma \sigma \mathrm{mn}}}{{J}_{\pi \pi \mathrm{mn}}+{J}_{\sigma \sigma \mathrm{mn}}}$$(20)
With the present numerical values we obtain P_{mn} ≈ 0.473 for all nearestneighbor orders. Even though the excitation wave is partially polarized, we obtain the same values of J_{σσmn} as above; both of the two mutually uncorrelated components of the incident field contribute to TEmode BSW excitation.
5.2 Sixway splitting
Let us consider a grating with hexagonal symmetry, which allows simultaneous excitation of six BSWs. The primitive vectors are now ${\mathit{a}}_{1}=d\widehat{\mathit{x}}$ and ${\mathit{a}}_{2}=(d/2)\widehat{\mathit{x}}+(\sqrt{3}d/2)\widehat{\mathit{y}}$, and the Wigner–Seitz primitive cell is a hexagon as shown in Figure 6a. It will, however, be convenient for our purposes to define a rectangular directlattice cell as in reference [20], which covers the spatial region −d/2 < x <d/2, $\sqrt{3}d/2<y<\sqrt{3}d/2$ in Figure 6a. This alternative lattice representation simplifies the visualization of the geometry. It also allows the use of FMM in a cartesian instead of a nonorthogonal basis, as in reference [20], though in the present work we actually used the latter basis.
Figure 6 (a) The spatial structure of a hexagonal grating. The hexagon shows the spatial Wigner–Seitz primitive cell, the green features illustrate the grating structure, and the blue rectangle shows the nonprimitive cartesian cell. (b) The spatialfrequency grid. The filled and empty dots represent the allowed and forbidden orders of the hexagonal lattice. The blue, yellow, and red circles have the same meaning as in Figure 2a. 
The spatialfrequency structure defined by the reciprocallattice primitive vectors ${\mathit{b}}_{1}=K\widehat{\mathit{x}}(K/\sqrt{3})\widehat{\mathit{y}}$, ${\mathit{b}}_{2}=(2K/\sqrt{3})\widehat{\mathit{y}}$ with K = 2π/d is illustrated in Figure 6b, where the solid green circles show the locations of the allowed orders in the cartesian (k_{x}, k_{y}) system. The empty circles represent the orders of the rectangular spatial lattice, which are forbidden by the hexagonal symmetry. The yellow circle connects the six nearest neighbors of the zeroth order that satisfy the condition for BSW excitation simultaneously: orders (m, n) = (+1, +1), (0, +2), (−1, +1), (−1, −1), (0, −2), (+1, −1) of the rectangular lattice, with exit planes at angles π/6 + qπ/3, q = 0,…,5. The excited BSWs propagate along the surface of the stack in these directions.
In hexagonal lattice geometry, the symmetry rules in equation (19) ensure J_{σσ,−1,−1} = J_{σσ,1,1}, J_{σσ,0,−2} = J_{σσ,0,2}, J_{σσ,−1,1} = J_{σσ,1,−1}. These symmetries leave us three pairs of orders to control, and we expect to need additional structural freedom compared to the 4wave case. Let us nevertheless see what designs are possible with circular pillars by following the same strategy as above. An important difference is that in the hexagonal geometry we do have structurally induced polarization conversion.
By optimizing r and h for pillars, we get h ≈ 0.422λ_{0} and r ≈ 0.254λ_{0}, which leaves a combined efficiency of ~0.884 available for the 6 orders of interest. The distribution of J_{σσmn} again depends on input polarization. We found that it is not possible to equalize the coupling exactly for all six orders, but using circularly polarized light with α = π/4 and δ = 0.486π we have J_{σσ,1,1} = 0.097, J_{σσ,0,2} = 0.102, and J_{σσ,1,−1} = 0.113, respectively. Similarly, for κ_{mn}, we have κ_{1,1} = 1.691, κ_{0,2} = 1.903, and κ_{1,−1} = 2.310. Though it is not of concern for the present purposes, it is worth noting that the diffraction efficiencies are: η_{1,1} ≈ 0.152, η_{0,2} ≈ 0.145, and η_{1,−1} ≈ 0.144. As with the fourwave case, the design with circular pillars works also for circularly polarized or unpolarized illumination but the exact values of J_{σσmn} depend on polarization, but are within the same range as above. For unpolarized illumination, the degrees of polarization of the individual orders are nearly equal, P_{mn} ≈ 0.3226.
In Figure 7, we show the field amplitude distribution associated with the sixway coupling geometry in the xzplane, i.e., crossing the multilayer, when illuminated with a 45° polarized light wave. The field is evaluated over 3unit cells, i.e., 3 grating periods, in the xdirection. It shows, as expected, a strong field on the upper medium. Such a structure is ideal for sensing applications, especially when providing multiple sensing areas thanks to the splitting of the BSW excitation.
Figure 7 Field amplitude distribution across the multilayer (xzplane). The illumination polarization was set to 45°. The dashed lines superimposed on the field represent the multilayer interfaces. 
6 Discussion
Throughout the paper we have considered normally incident illumination. The use of nonnormal incidence could potentially allow us to consider other combinations of diffracted orders being simultaneously resonant. Changing the angle of incidence moves the grid of diffracted orders transversely in Figure 2a with respect to the circles centered at the origin. For instance, if the propagation direction of the incident field is chosen as (k_{xi}, k_{yi}) = (0, k_{yi}), increasing k_{yi} moves the grid downwards in k_{y} direction, giving the orders at positions k_{xm} = mK_{x}, k_{yn} = k_{yi} +nK_{y}. Hence the three orders (m, n) = (0, 0) and (m, n) = (±1, 0) would have a common radial spatial frequency if ${k}_{y\mathrm{i}}=({K}_{x}^{2}+{K}_{y}^{2})/2{K}_{y}$, being therefore available for 3way BSW excitation. To avoid order (0, 2) from occupying the same ring as the zeroth order, we would need to choose K_{x} ≠ K_{y}. However, placing the (yellow) BSW resonance ring outside the blue ring in Figure 2a requires k_{ρ00} > k_{0}n_{sub}, which is not possible with incidence from air. Hence a Kretschmann excitation geometry would be needed, thus sacrificing the compact footprint of the setup.
As an alternative to the geometry considered in this paper, we could consider having the splitter grating and the BSW stack as an integrated structure. This would still allow a compact platform at normal incidence, but the grating design and BSW stack design would not be independent anymore. As a first drawback, the splitter grating would most likely have to be rather thick (~λ_{0}) to suppress the zeroth transmitted order, preventing the possibility of etching it in the top ML layer. As a consequence, a strong degrading effect in the excited BSWs would be expected. Alternatively, one could use a highly indexmodulated splitter grating with a flat top surface immediate below the stack. This would partially alleviate the dependency in the ML and the grating design, but presumably the BSWs would be less affected.
7 Conclusions
In summary, we considered grating design for 2way, 4way, and 6way BSW coupling at normally incident but arbitrarily polarized illumination of gratings with linear, square, and hexagonal symmetries. The planewave designs feature ideal TEmode BSW coupling in the twowave case. In the other cases the nonresonant parts of the excitation orders cannot be eliminated simultaneously, and they actually dominate the resonant (σ polarized) parts by a factor of ~2.8 in the 4wave case. The 6wave case reveals the opposite observation with the resonant part dominating by a factor of ~1.968 making them ideal candidate for TEmode BSW excitation. Nevertheless, in the 4wave case all four coupling ratios can be made equal, and in the 6wave case practically equal, for several input polarization states of practical significance.
Acknowledgments
We want to pay our respects to our mentor, teacher, and friend, Dr Jari Turunen, who passed away in May 2023. He was a dedicated professor in the Center for Photonics Sciences of the Department of Physics and Mathematics at the University of Eastern Finland, Joensuu, who was passionate about research.
Funding
The work was partially funded by the Academy of Finland through project 333938 and the Flagship Programme PREIN (346518). E. Descrovi acknowledges the funding received by Italian “Ministero dell’Universitá e della Ricerca” under the “Dipartimento di Eccellenza 20182022” program.
Conflicts of Interest
The authors declare no conflicts of interest.
Data availability statement
All data generated or analyzed during this study are included in this published article.
Author contribution statement
ALA developed the codes and performed the numerical simulations. HP contributed in developing the codes and to the discussion. ED and MR supervised the work. MR and JT proposed the concept. JT developed the theory. All authors contributed to the manuscript preparation, have read and approved the final manuscript.
References
 Robertson W.M., May M.S. (1999) Surface electromagnetic wave excitation on onedimensional photonic bandgap arrays, Appl. Phys. Lett. 74, 1800–1802. [CrossRef] [Google Scholar]
 Liscidini M., Sipe J.E. (2007) Enhancement of diffraction for biosensing applications via Bloch surface waves, Appl. Phys. Lett. 91, 253125. [NASA ADS] [CrossRef] [Google Scholar]
 Villa F., GasparArmenta J.A. (2004) Photonic crystal to photonic crystal surface modes: narrowbandpass filters, Opt. Express 12, 2338–2355. [NASA ADS] [CrossRef] [Google Scholar]
 Mogni E., Pellegrini G., GilRostra J., Yubero F., Simone G., Fossati S., Dostälek J., Martínez Vázquez R., Osellame R., Celebrano M., Finazzi M., Biagioni P. (2022) Onedimensional photonic crystal for surface mode polarization control, Adv. Opt. Mater. 10, 21, 2200759. [CrossRef] [Google Scholar]
 Asilevi A.L., Pesonen H., Pelisset S., Descrovi E., Roussey M., Turunen J. (2022) Pulse modulation by Bloch surface wave excitation, Opt. Lett 47, 2574–2577. [NASA ADS] [CrossRef] [Google Scholar]
 Descrovi E., Frascella F., Sciacca B., Geobaldo F., Dominici L., Michelotti F. (2007) Coupling of surface waves in highly defined onedimensional porous silicon photonic crystals for gas sensing applications, Appl. Phys. Lett. 91, 241109. [CrossRef] [Google Scholar]
 Gulkin D.N., Popkova A.A., Afinogenov B.I., Shilkin D.A., Kuršelis K., Chichkov B.N., Bessonov V.O., Fedyanin A.A. (2021) Miedriven directional nanocoupler for Bloch surface wave photonic platform, Nanophotonics 10, 2939–2947. [Google Scholar]
 Wang M., Zhang H., Kovalevich T., Salut R., Kim M.S., Suarez M.A., Bernal M.P., Herzig H.P., Lu H., Grosjean T. (2018) Magnetic spin–orbit interaction of light, Light Sci. Appl. 7, 24. [NASA ADS] [CrossRef] [Google Scholar]
 Wang R., Xia H., Zhang D., Chen J., Zhu L., Wang Y., Yang E., Zang T., Wen X., Zou G., Wang P. (2017) Bloch surface waves confined in one dimension with a single polymeric nanofibre, Nat. Commun. 8, 14330. [Google Scholar]
 Safronov K.R., Bessonov V.O., Akhremenkov D.V., Sirotin M.A., Romodina M.N., Lyubin E.V., Soboleva I.V., Fedyanin A.A. (2022) Miniature Otto prism coupler for integrated photonics, Laser Photonics Rev. 16, 2100542. [NASA ADS] [CrossRef] [Google Scholar]
 Bezus E.A., Bykov D.A., Doskolovich L.L. (2021) Integrated diffraction gratings on the Bloch surface wave platform supporting bound states in the continuum, Nanophotonics 10, 4331–4340. [CrossRef] [Google Scholar]
 Kovalevich T., Boyer P., Suarez M., Salut R., Kim M.S., Herzig H.P., Bernal M.P., Grosjean T. (2017) Polarization controlled directional propagation of Bloch surface wave, Opt. Express 25, 5710. [NASA ADS] [CrossRef] [Google Scholar]
 Scaravilli M., Micco A., Castaldi G., Coppola G., Gioffrë M., Iodice M., La Ferrara V., Galdi V., Cusano A. (2018) Excitation of Bloch surface waves on an optical fiber tip, Adv. Opt. Mater. 6, 1800477. [CrossRef] [Google Scholar]
 Koju V., Robertson W.M. (2017) Leaky Blochlike surface waves in the radiationcontinuum for sensitivity enhanced biosensors via azimuthal interrogation, Sci. Rep. 7, 3233. [NASA ADS] [CrossRef] [Google Scholar]
 Scaravilli M., Castaldi G., Cusano A., Galdi V. (2016) Gratingcouplingbased excitation of Bloch surface waves for labonfiber optrodes, Opt. Express 24, 27771. [CrossRef] [Google Scholar]
 Wang R., Lei X., Jin Y., Wen X., Du L., Wu A., Zayats A.V., Yuan X. (2021) Directional imbalance of Bloch surface waves for ultrasensitive displacement metrology, Nanoscale 13, 11041–11050. [Google Scholar]
 Deng C.Z., Ho Y.L., Yamahara H., Tabata H., Delaunay J.J. (2022) Nearzeroindex slabs on bloch surface wave platform for longrange directional couplers and optical logic gates, ACS Nano 16, 2224–2232. [CrossRef] [Google Scholar]
 Lei X., Wang R., Liu L., Xu C., Wu A., Zhan Q. (2022) Multifunctional onchip directional coupler for spectral and polarimetric routing of Bloch surface wave, Nanophotonics 11, 21, 4627–4636. [Google Scholar]
 Kim H., Park J., Lee B. (2012) Fourier Modal Method and it Applications in Computational Nanophotonics, CRC Press. [Google Scholar]
 Noponen E., Turunen J. (1994) Eigenmode method for electromagnetic synthesis of diffractive elements with threedimensional profiles, J. Opt. Soc. Am. A 11, 2494–2502. [NASA ADS] [CrossRef] [Google Scholar]
 Safronov K.R., Gulkin D.N., Antropov I.M., Abrashitova K.A., Bessonov V.O., Fedyanin A.A. (2020) Multimode interference of Bloch surface electromagnetic waves, ACS Nano 14, 10428–10437. [CrossRef] [Google Scholar]
 Dubey R., Lahijani B.V., Barakat E., Häyrinen M., Roussey M., Kuittinen M., Herzig H.P. (2016) Nearfield characterization of a Blochsurfacewavebased 2D disk resonator, Opt. Lett. 41, 4867–4870. [NASA ADS] [CrossRef] [Google Scholar]
 Tervo J., Turunen I.A., Bai B. (2008) A general approach to the analysis and description of partially polarized light in rigorous grating theory, J. Eur. Opt. Soc: Rapid Publ. 3, 08004. [NASA ADS] [CrossRef] [Google Scholar]
 Mandel L., Wolf E. (1995) Optical Coherence and Quantum Optics, Cambridge. [CrossRef] [Google Scholar]
All Figures
Figure 1 Principle of MBSW generation: the twobeam case. A binary linear surfacerelief grating defined by period d, ridge width c, and ridge height h on the bottom of a substrate of thickness H splits the input beam into two diffracted orders m = −1 and m = +1, which excite BSWs on the top surface of the substrate by interaction with the multilayer stack (ML). We assumed n_{sub} = 1.462 and n_{sup} is air. 

In the text 
Figure 2 (a) Diffraction orders of a rectangular lattice in spatialfrequency representation at normal incidence. Diffraction orders are represented by dots at positions k_{xm} = mK_{x}, k_{y} = nK_{y}. Blue and red circles represent the cutoff radial spatial frequencies k_{ρ} = k_{0}n_{sub} and k_{ρ} = k_{0}n_{sup}, respectively, between which BSW excitation is possible. The yellow circle indicates the radial spatial frequency of BSW for a given ML. (b) Definition of the propagation angles (θ_{mn}, ϕ_{mn}) of a single transmitted diffracted order (m, n) in its exit plane (the grey rectangle) and the π − σ basis of the diffracted electric field. 

In the text 
Figure 3 (a) Definition of the multilayer structure and notation. (b) Stack parameters as a function of the ratio k_{x}/k_{0} = n_{eff} for TEmode BSW excitation with N = 6 bilayers: h_{H}/λ_{0} (blue), h_{L}/λ_{0} (red), and h_{T}/λ_{0} (black). The dots mark the position k_{x}/k_{0} = 1.1209 for BSW excitation at 50° angle of incidence. 

In the text 
Figure 4 Design of twoway beam splitters. (a) Optimum values of the fill factor f (red) and the relief depth h/λ_{0} (blue) as a function of the exit angle of the first diffracted order in TE polarization. (b) The corresponding firstorder diffraction efficiency η_{+1,0} = η_{−1,0} (black), the efficiency η_{0,0} of the zeroth transmitted order (red), and that of the zeroth reflected order (blue), which is 1 − 2η_{+1,0} − η_{0,0} due to energy conservation. The inset shows the grating structure and direction of illumination. 

In the text 
Figure 5 Response of the full structure (grating, substrate, multilayer and superstrate). (a) and (b) Reflected (black curves) and transmitted (red curves) first diffracted orders as a function of the normalized wavevector (a) and wavelength (b). 

In the text 
Figure 6 (a) The spatial structure of a hexagonal grating. The hexagon shows the spatial Wigner–Seitz primitive cell, the green features illustrate the grating structure, and the blue rectangle shows the nonprimitive cartesian cell. (b) The spatialfrequency grid. The filled and empty dots represent the allowed and forbidden orders of the hexagonal lattice. The blue, yellow, and red circles have the same meaning as in Figure 2a. 

In the text 
Figure 7 Field amplitude distribution across the multilayer (xzplane). The illumination polarization was set to 45°. The dashed lines superimposed on the field represent the multilayer interfaces. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.