Shortcuts to adiabaticity in waveguide couplers–theory and implementation

ABSTRACT

Adiabatic processes are ubiquitous in physics and engineering. A drawback of such processes is that they tend to be slow, either in time–for atomic systems, for example–or in space–for photonic systems. A number of techniques have been developed over the years, generically referred to as Shortcuts to Adiabaticity, that promise to speed up the evolution of the process without compromising performance. Here we review and compare these techniques, and evaluate their performance using full numerical simulations of realistic two-waveguide couplers, which perform a key function in photonic circuits.

Graphical Abstract

KEYWORDS:

• Adiabatic devices
• waveguide couplers

1. Introduction

Adiabatic processes, first considered by Fock and Born [1], arise in areas of physics and engineering such as atomic physics, quantum mechanics and photonics. Adiabatic processes occur slowly enough that the wave function or field remains in the same Eigenstate or mode throughout the evolution, even when the physical nature of this Eigenstate changes strongly [2–6]. The process may evolve as a function of time, as usually applies in quantum mechanics, or of propagation length, as typically applies in photonics. Here we are interested in photonics and in particular adiabatic waveguide couplers, which aim to couple light as efficiently as possible from one waveguide to another. Since the formalism for quantum mechanical and photonics processes is identical, we will take the evolution to be in $z$. Thus, when we refer to the length of a device, in an atomic physics context it refers to the duration of the process.

Most photonic systems require the efficient transfer of light from one waveguide into another [7], so a number of different approaches have been developed to achieve it: one of these is butt coupling, in which the two waveguides are placed end-to-end [8]. Such devices are by their nature very short, but the efficiency of butt-coupling can be low when the modal profiles and the propagation constants of the modes in the two waveguides are poorly matched. Another approach is the directional coupler, illustrated in Fig. 1(a), in which the two waveguides are placed side-by-side. If the two individual waveguides have the same propagation constants and the light is coupled into an equal superposition of the two supermodes of the coupled waveguides, the light couples back and forth periodically over the length of the device (Fig. 1(b)). If the device is terminated at the correct position, all light exits through the other waveguide. Though directional couplers can be fairly short, the requirement that the waveguides have the same propagation constants, and the dependence of the coupling length on wavelength, constrain the performance in terms of fabrication tolerances and bandwidth [9,10].

In contrast to these two devices, adiabatic couplers, introduced into photonics by Cook [11], are relatively robust devices that can have a large bandwidth [6,9,12–18]. The type of adiabatic coupling device we are considering here is illustrated in Fig 1(c). We consider the two lowest supermodes of the device; approximately each supermode can be thought of as a superposition of the fundamental modes of the isolated waveguides 1 and 2. However, the device is designed such that at $z=0$ and $z=\mathrm{\ell }$ the supermodes closely match the modes of the individual waveguides: for the device in Fig. 1(c), the fundamental supermode at $z=0$ corresponds to the isolated mode of waveguide 1, whereas at $z=\mathrm{\ell }$ it corresponds to that of waveguide 2. For the device in the schematic, this is achieved by varying the widths of the waveguides. In a perfectly adiabatic process, light that is coupled into the fundamental supermode enters through waveguide 1, remains in the fundamental supermode, and exits via waveguide 2, ideally with 100% efficiency. A consequence of adiabaticity is that, provided that the supermodes coincide with each of the individual waveguides at $z=0$ and $z=\mathrm{\ell }$, the device performance is relatively insensitive to the perturbations of the device or even to the details of the device design [9,17].

Figure 1. (a) Schematic of a directional coupler: two supermodes, with different propagation constants, are excited at the input in waveguide 1. (b) The supermodes’ interference leads to a periodically varying output with length $\mathrm{\ell }$. (c) Schematic of an adiabatic coupler: Light is input in waveguide 1, corresponding to a supermode at $z=0$. The light remains in this mode, which at the end of the device corresponds to waveguide 2. The dashed line schematically indicates the propagation of the light. (d) As the device length $\mathrm{\ell }$ increases, the crosstalk $\mu$ has discrete zeros under an envelope that decreases with $\mathrm{\ell }$

Even though adiabatic devices promise excellent performance, because of the requirement of adiabaticity they tend to be quite long–ideally they are infinitely long! It is therefore not surprising that a number of techniques have been developed that promise to minimize this length, while maintaining performance. Most of these techniques are referred to generically as Shortcuts to Adiabaticity (STA) [3–6]—developed initially in the context of quantum mechanics, these have been considered for applications in optics since 2009 [19]. In Section 3, we discuss these techniques and clarify their mutual relationships.

One shortcut technique we discuss in Section 3 is based on the Lewis-Riesenfeld invariant theory [20–25] in which the device design is reverse-engineered to accelerate adiabatic processes, and it has been used in applications such as atomic transport [26–28] and trap compression and expansions [29–32]. Other STA methods include the counter-diabatic (CD) approach developed by Demirplak and Rice [33] and the related transitionless tracking method of Berry [34]. Starting from an initial device, they introduce an additional contribution that cancels the diabatic coupling in the original device, reducing the effect of the imperfect adiabaticity. This approach has been applied to speed up rapid adiabatic passage (RAP) in two level atoms [35] and was experimentally implemented in Bose-Einstein condensates in optical lattices [36]. Alternatively, in the fast quasi-adiabaticity (FAQUAD) [37] technique, the adiabaticity is violated to the same extent throughout the device. We also discuss the relationships between these methods, first addressed by Chen et al [21]. Although we have included a large set of references, Figure 2 summarises some of the main references–the foundational papers, major reviews, key references applying STA techniques in photonics–and provides a timeline.

Figure 2. Timeline with foundational papers, major reviews, and key references applying STA techniques in photonics. The figure is not intended to be exhaustive. Key references have a coloured background

While developed in the context of quantum mechanics and atomic physics, these techniques can also be applied to photonic devices when these are described using coupled mode theory (CMT). This theory enables a general, simple yet powerful description of photonic devices in which a complex physical design is expressed in terms of a reduced number of parameters. Using this framework, a coupler such as in Figure 1 can be described in terms of only two real functions, corresponding to the difference in propagation constant of the modes of the device ($\mathrm{\Delta }$), and the strength of the coupling between the waveguides ($\kappa$), each of which depends on longitudinal position, and thus $\mathrm{\Delta }=\mathrm{\Delta }(z)$ and $\kappa =\kappa (z)$.

The aim of this Review is to examine the adiabatic literature in the context of waveguide couplers, and to discuss which of the techniques mentioned above can be useful in this context. Although the techniques have been applied to other devices such as junctions [38,39], photonic lattices [40], demultiplexers [41], mode converters [42] and filters [43], curved waveguides [19], as well as multimode waveguides [44] and interference splitters [45], we choose to consider the coupling between two waveguides [46–49] (as in Figure 1) as this is a simple nontrivial case with only two relevant modes. Therefore, the Hamiltonian that describes the interaction between these modes can be written as a $2\times 2$ matrix with position-dependent coefficients.

We consider lossless systems and we take this $2\times 2$ matrix to be real and symmetric, so that it has real Eigenvalues and orthogonal Eigenvectors. This choice deserves two comments. The first indicates a difference with quantum mechanics, where the $2\times 2$ matrix is generally Hermitian, but in the description of the coupling between waveguides the matrix cannot have complex elements, and thus the description by a Hermitian matrix reduces to that of a real, symmetric one. The second point is that for strongly coupled waveguides, the interaction between lossless waveguides requires a real matrix that is asymmetric [50,51]. Such matrices are of course not Hermitian, and cannot therefore be described by the theory we are reviewing here, though we note that STA has been generalized to non-Hermitian systems [52,53].

The outline of this Review is as follows. In Section 2 we give a general overview of adiabatic couplers and their mathematical description. In Section 3 we review the key shortcut methods that have been most widely used and discuss their mutual relationships. In Section 4 we illustrate the use of these methods in a physically realistic geometry and compare to alternative adiabatic couplers generated by other techniques. Then, in Section 5 we compare the sensitivity of a number of different adiabatic couplers to the presence of noise in the device parameters, motivated to evaluate robustness to fabrication errors. Finally, in Section 6 we discuss our findings and conclude.

2. Adiabatic processes and FAQUAD

Waveguide couplers as in Figure 1 must be described by several parameters (refractive indices and dimensions of each waveguide, spacing of the waveguides, refractive index of the region between and outside the waveguides), all of which are functions of propagation length. As discussed in Section 1, in coupled mode theory only two parameters are required. In coupled mode theory the field in a waveguide is written as the product of a mode field, which depends on the transverse coordinates, and a $z$-dependent factor describing the propagation. For an isolated waveguide this $z$-dependent factor is $F\mathrm{\infty }exp(i\tilde{\beta }z)$, where $F$ is the mode amplitude, defined such that $|F{|}^2$ corresponds to the power carried by the mode. and $\tilde{\beta }$ is the wave number of the mode. The propagation in the isolated waveguides is thus described by $idF/dz=-\tilde{\beta }F$.

If two waveguides $1$ and $2$ are brought together, and allowed to interact linearly then we find

(1) $\begin{array}{c}i\frac{d{F}_1}{dz}=-{\beta }_1{F}_1+\kappa F_2,\\ i\frac{d{F}_2}{dz}=+\kappa F_1-{\beta }_2{F}_2,\end{array}$(1)

where $\kappa$ represents the coupling between the waveguides and both $\kappa$ and the ${\beta }_{1,2}$ are real functions of $z$. We note that ${\tilde{\beta }}_{1,2}\ne {\beta }_{1,2}$ due to the proximity of the other waveguide. Now defining ${\beta }_{1,2}=\bar{\beta }\mp \mathrm{\Delta }$ and $G=Fexp(-i\bar{\beta }z)$, Equation (1)(1) $\begin{array}{c}i\frac{d{F}_1}{dz}=-{\beta }_1{F}_1+\kappa F_2,\\ i\frac{d{F}_2}{dz}=+\kappa F_1-{\beta }_2{F}_2,\end{array}$(1) can be rewritten as

(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2)

which defines the ‘Hamiltonian’ $H$, and where the ket $|G⟩$ has the elements $G_{1,2}$.

In the following analysis of these equations we follow Louisell [46], who, to the best of our knowledge, reported the first systematic investigation of adiabatic processes in the context of photonics. The instantaneous Eigenvalues of $H$ are $\pm \mathrm{\Gamma }\equiv \pm \sqrt{{\mathrm{\Delta }}^2+{\kappa }^2}$ with associated Eigenvectors

(3) $|+⟩=\left(\begin{array}{c}cos\theta /2\\ sin\theta /2\end{array}\right),\mathrm{a}\mathrm{n}\mathrm{d}|-⟩=\left(\begin{array}{c}sin\theta /2\\ -cos\theta /2\end{array}\right),$(3)

corresponding to the instantaneous supermodes of the device, and where

(4) $tan\theta =\kappa /\mathrm{\Delta }.$(4)

We now write the general solution to Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) as $|G⟩=w_{+}|+⟩+w_{-}|-⟩$; when $\mathrm{\Delta }$ and $\kappa$ are constant, then the $w_{\pm }$ vary harmonically. In general though, $\mathrm{\Delta }$ and $\kappa$ depend on $z$, and the $w_{\pm }$ satisfy the coupled equations [46]

(5) $\begin{array}{c}\frac{d{w}_{+}}{dz}+i\mathrm{\Gamma }{w}_{+}=-\frac{1}{2}\frac{d\theta }{dz}{w}_{-},\\ \frac{d{w}_{-}}{dz}-i\mathrm{\Gamma }{w}_{-}=+\frac{1}{2}\frac{d\theta }{dz}{w}_{+}.\end{array}$(5)

Equations (4)(4) $tan\theta =\kappa /\mathrm{\Delta }.$(4) and (5(5) $\begin{array}{c}\frac{d{w}_{+}}{dz}+i\mathrm{\Gamma }{w}_{+}=-\frac{1}{2}\frac{d\theta }{dz}{w}_{-},\\ \frac{d{w}_{-}}{dz}-i\mathrm{\Gamma }{w}_{-}=+\frac{1}{2}\frac{d\theta }{dz}{w}_{+}.\end{array}$(5) ) confirm that when $\mathrm{\Delta }$ and $\kappa$ are constant, the $w_{\pm }$ vary as $e^{\pm i\mathrm{\Gamma }z}$ and do not couple. Variations in $\mathrm{\Delta }$ and $\kappa$ that cause variations in $\theta$, lead to coupling of the $w_{\pm }$. In the truly adiabatic limit the derivatives on the right-hand sides in Equations (5)(5) $\begin{array}{c}\frac{d{w}_{+}}{dz}+i\mathrm{\Gamma }{w}_{+}=-\frac{1}{2}\frac{d\theta }{dz}{w}_{-},\\ \frac{d{w}_{-}}{dz}-i\mathrm{\Gamma }{w}_{-}=+\frac{1}{2}\frac{d\theta }{dz}{w}_{+}.\end{array}$(5) are arbitrarily small and the equations again decouple.

Adiabatic coupling can then be understood by taking $\theta (0)=0$ and $\theta (\mathrm{\ell })=\pi$. If light is coupled into waveguide 1, say, at $z=0$, it is in supermode $|+⟩$. In an adiabatic coupler without supermode coupling, the light remains in this supermode, and therefore exits entirely through waveguide 2. Of course, for the light to remain rigorously in the same supermode throughout its propagation, the device needs to be infinitely long. All the techniques discussed in this Review correspond to strategies to minimise the device length without compromising device performance.

Although the parameters that enter the Hamiltonian are $\mathrm{\Delta }$ and $\kappa$, they enter as the combinations $\mathrm{\Gamma }$ and $\theta$ in the solutions. However, only variations in $\theta$ cause changes in the supermodes; thus, if $\mathrm{\Delta }$ and $\kappa$ vary at the same rate, so $\theta$ is constant, then the supermodes are unchanged.

Equations (5)(5) $\begin{array}{c}\frac{d{w}_{+}}{dz}+i\mathrm{\Gamma }{w}_{+}=-\frac{1}{2}\frac{d\theta }{dz}{w}_{-},\\ \frac{d{w}_{-}}{dz}-i\mathrm{\Gamma }{w}_{-}=+\frac{1}{2}\frac{d\theta }{dz}{w}_{+}.\end{array}$(5) can be solved formally at various levels of approximation. Assuming that the field amplitude in the original super mode remains close to unity, the crosstalk $\mu$, that is, the fraction of the incoming power that is not coupled over at $z=\mathrm{\ell }$ can be written as [46]

(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6)

where $d\rho /dz=\mathrm{\Gamma }$, and $\rho$ is an effective position. Note from Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) that $\mu$ is, in essence, the Fourier transform of $d\theta /dz$, as is particularly clear when $\mathrm{\Gamma }$ is constant. Now recall the theorem that the Fourier transform of a function with a discontinuous $m^{\mathrm{t}\mathrm{h}}$ derivative, asymptotically is a power law with exponent $-(m+1)$ [47,54]. Applying this to Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) we conclude that the asymptotic behaviour of $\mu$, i.e. the behaviour as $\mathrm{\ell }\to \mathrm{\infty }$, is determined by the inevitable discontinuities of $\theta$ or any of its derivatives at $z=0$ and $z=\mathrm{\ell }$ (assuming that $\theta$ and all of its derivatives are smooth elsewhere): when the $m^{\mathrm{t}\mathrm{h}}$ derivative is discontinuous, $\mu \propto {\mathrm{\ell }}^{-2(m+1)}$ [47].

Figure 1(d) shows a typical result for the crosstalk versus device length as calculated by integrating the coupled mode Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) . Note that the crosstalk has discrete zeros, superimposed on a monotonically decreasing envelope. Before discussing the results in Figure 1(d), we note that Equations (5)(5) $\begin{array}{c}\frac{d{w}_{+}}{dz}+i\mathrm{\Gamma }{w}_{+}=-\frac{1}{2}\frac{d\theta }{dz}{w}_{-},\\ \frac{d{w}_{-}}{dz}-i\mathrm{\Gamma }{w}_{-}=+\frac{1}{2}\frac{d\theta }{dz}{w}_{+}.\end{array}$(5) and (6(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) ) show that the relevant quantity that determines the degree of adiabaticity involves $d\theta /dz$; in fact, the adiabaticity criterion can be written as [46]

(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7)

Given inequality (7), one strategy to designing a coupler is to take $\eta$ to be constant–this is the fast quasi-adiabatic (FAQUAD) approach, which was introduced formally by Martínez-Garaot et al. [37] and applied by Hung et al. [55], but was used prior to that by Sun et al. [56]. From Equations (5)(5) $\begin{array}{c}\frac{d{w}_{+}}{dz}+i\mathrm{\Gamma }{w}_{+}=-\frac{1}{2}\frac{d\theta }{dz}{w}_{-},\\ \frac{d{w}_{-}}{dz}-i\mathrm{\Gamma }{w}_{-}=+\frac{1}{2}\frac{d\theta }{dz}{w}_{+}.\end{array}$(5) , (6(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) ) and (7(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) ), decreasing $\eta$ leads to lower crosstalk, but with a longer device. Hence, the choice is a compromise between performance and device length. Note that within the constraint of Equation (7)(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) , any $\theta (z)$ that is continuous and that has a continuous first derivative lends itself to the FAQUAD approach. Although not obvious from Equation (7)(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) , the FAQUAD procedure has the advantage that it can be applied without the need to use a coupled mode description of the device.

To illustrate the results in this section, we consider a set of couplers with parameters that can be varied systematically [47]: we consider the elements of the Hamiltonian to take the values

(8) $\begin{array}{rl}& \mathrm{\Delta }=\mathrm{\Gamma }cos(\pi P_n(\zeta )),\\ & \kappa =\mathrm{\Gamma }sin(\pi P_n(\zeta )),\end{array}$(8)

where the normalised propagation length $\zeta =z/\mathrm{\ell }$, with $\mathrm{\ell }$ the device length, and thus $0\le \zeta \le 1$, and $\mathrm{\Gamma }$ is constant. The $P_n(\zeta )$ are polynomials which satisfy $P_n(\zeta )+P_n(1-\zeta )=1$, with $P_n(0)=0$ and thus $P_n(1)=1$. They are most easily defined through their derivative (with respect to the argument) through [47]

(9) $P{{ }^{\mathrm{\prime }}}_n(\zeta )=\frac{(2n+1)!}{{(n!)}^2}(\zeta (1-\zeta ){)}^n,$(9)

and thus $P_0(\zeta )=\zeta$ and $P_1(\zeta )=3{\zeta }^2-2{\zeta }^3$, etc. As illustrated in Figure 3(a), at $\zeta =0,1$, the lowest $n$ derivatives of these polynomials are continuous, whereas higher ones are discontinuous. All of the designs based on Equation (8)(8) $\begin{array}{rl}& \mathrm{\Delta }=\mathrm{\Gamma }cos(\pi P_n(\zeta )),\\ & \kappa =\mathrm{\Gamma }sin(\pi P_n(\zeta )),\end{array}$(8) have the property that $\theta$ increases smoothly from $\theta =0$ at $\zeta =0$ to $\theta =\pi$ at $\zeta =1$. As $n$ increases, the discontinuities at the edges become weaker, and so $\theta$ needs to increase more rapidly in the central part of the coupler. Weakening the discontinuities at the edges improves the performance for large $\mathrm{\ell }$.

Figure 3. Three different polynomial couplers $n=-2,0,2$ at different stages of the design. (a) shows the mixing angle $\theta$ along the normalised length of the device. $P_{-\mathrm{\infty }}$ (a directional coupler) and $P_{\mathrm{\infty }}$ are also plotted for comparison. (b) Normalised Hamiltonian components, with coupling coefficient $\kappa$ (solid) and phase velocity mismatch $\mathrm{\Delta }$ (dashed) versus normalised position along the device. (c) Crosstalk from solving Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) (solid lines) for different device lengths, keeping $\mathrm{\Gamma }$ constant. The crosstalk envelope (dashed line) decreases for increasing $n$

To see the drawback of increasing $n$ we consider the expression for $\mu$ in more detail. For constant $\mathrm{\Gamma }$ and for the symmetric devices we are considering here, the crosstalk vanishes when

(10) ${\int }_0^{\frac{1}{2}}\frac{d\bar{\theta }}{d\zeta }cos(2\mathrm{\Gamma }\mathrm{\ell }\zeta )d\zeta =0,$(10)

where $\bar{\theta }(\zeta )=\theta (\zeta +1/2)$. For large $n$, $d\theta /d\zeta$ is concentrated near $\zeta =1/2$, and so the first zero occurs for large values of $\mathrm{\Gamma }\mathrm{\ell }$, whereas for small $n$ the first zero occurs at smaller values of $\mathrm{\Gamma }\mathrm{\ell }$. To illustrate this, for $n=0$, the first zero in the crosstalk occurs when $\mathrm{\Gamma }\mathrm{\ell }=\pi$, whereas for $n=1$ it occurs at $\mathrm{\Gamma }\mathrm{\ell }=1.43\pi$, and for $n=2$ at $\mathrm{\Gamma }\mathrm{\ell }=1.83\pi$. By comparison, for a directional coupler, for which $\mathrm{\Delta }=0$ and thus $\kappa =\mathrm{\Gamma }$, the first zero occurs when $\mathrm{\Gamma }\mathrm{\ell }=\pi /2$.

The couplers defined by Equations (8)(8) $\begin{array}{rl}& \mathrm{\Delta }=\mathrm{\Gamma }cos(\pi P_n(\zeta )),\\ & \kappa =\mathrm{\Gamma }sin(\pi P_n(\zeta )),\end{array}$(8) and (9(9) $P{{ }^{\mathrm{\prime }}}_n(\zeta )=\frac{(2n+1)!}{{(n!)}^2}(\zeta (1-\zeta ){)}^n,$(9) ) thus illustrate the practical difficulty in choosing an ‘optimal’ design. For the polynomials with large $n$, the envelope (see Figure 3(c)) decreases rapidly with $\mathrm{\Gamma }\mathrm{\ell }$ [47]. Note that the crosstalk envelope represents the degree to which the mode remains in its Eigenstate during propagation. It thus indicates robustness to variations in the ideal design. For larger $n$, the first zero in the crosstalk appears at longer lengths. In contrast, for small $n$ the envelope decreases slowly, but the first zero appears at shorter lengths. However, in that case the design competes with directional couplers. Thus, the criterion that is used to optimize the coupler has a strong influence on its eventual shape.

For the coupler with $n=0$, the parameter $\eta$ defined in Equation (7)(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) is constant, and it thus satisfies the FAQUAD condition. To compare its performance with other couplers, we extend the polynomial Ansatz from Equations (8)(8) $\begin{array}{rl}& \mathrm{\Delta }=\mathrm{\Gamma }cos(\pi P_n(\zeta )),\\ & \kappa =\mathrm{\Gamma }sin(\pi P_n(\zeta )),\end{array}$(8) and (9(9) $P{{ }^{\mathrm{\prime }}}_n(\zeta )=\frac{(2n+1)!}{{(n!)}^2}(\zeta (1-\zeta ){)}^n,$(9) ) to negative integers $n$. We do so by defining $P_{-n}(\zeta )=P_n^{-1}(\zeta )$, where the superscript indicates the inverse of the function (rather than its reciprocal). Thus, for all integers $n$, the properties listed between Equations (8)(8) $\begin{array}{rl}& \mathrm{\Delta }=\mathrm{\Gamma }cos(\pi P_n(\zeta )),\\ & \kappa =\mathrm{\Gamma }sin(\pi P_n(\zeta )),\end{array}$(8) and (9(9) $P{{ }^{\mathrm{\prime }}}_n(\zeta )=\frac{(2n+1)!}{{(n!)}^2}(\zeta (1-\zeta ){)}^n,$(9) ) are maintained. In this way, we obtain a set of couplers that range from those that are highly discontinuous at the edges, and slowly varying in the centre ($n$ large and negative), to couplers which are smooth at the edges and rapidly varying in the centre ($n$ large and positive), with the FAQUAD coupler as a compromise between these extremes. We discuss these further in Sections 4 and 5.

3. Shortcuts to adiabaticity: methods and relations

In this section, we describe the three STA methods that have found most widespread use [57]: the invariant method of Lewis and Riesenfeld [20], the work of Demirplak and Rice in counter-diabatic (CD) driving [33], and Berry’s transitionless tracking (TT) algorithm [34]. We first give general descriptions of each of these methods at the operator level and then explore how they can be applied to solve the particular problem of power transfer in a waveguide coupler. Furthermore, we discuss connections between these methods. To better facilitate this, notation is kept consistent between subsections where applicable. Overall, these methods provide a way to start and end in a pure state of the Hamiltonian for a given device length. We explicitly show the $z$-dependence in all equations in this section for clarity.

3.1. Lewis-Riesenfeld invariants

The original work of Lewis and Riesenfeld [20] proposes to seek a Hermitian invariant $I(z)$ for a dynamic system driven by $H(z)$ satisfying

(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11)

for all $z$. Assuming the basis set is complete since the invariant arises from some observables, we can write

(12) $I(z)=\sum _n\left|{\varphi }_n(z)\right.⟩{\lambda }_n\left.⟨{\varphi }_n(z)\right|,$(12)

with $\left|{\varphi }_n(z)\right.⟩$ the orthonormal Eigenstates and ${\lambda }_n$ the corresponding Eigenvalues. By considering the $z$-derivative of this expression in conjunction with Equation (11)(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11) , Lewis and Riesenfeld [20] show that ${\lambda }_n$ must be constant. In general, the method uses ‘reverse engineering,’ with an assumed invariant, and then derives constraints on a Hamiltonian such that Equation (11)(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11) is satisfied. In doing so, states can be represented by

(13) $|G(z)⟩=\sum _n{c}_n{e}^{i{\alpha }_n(z)}\left|{\varphi }_n(z)\right.⟩,$(13)

where the $c_n$ are $z$-independent amplitudes and ${\alpha }_n(z)$ is the Lewis-Riesenfeld phase for state $n$. By requiring that the right hand side satisfies the dynamics in Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) for any choice of amplitudes, Lewis and Riesenfeld [20] show that

(14) $\frac{d{\alpha }_n}{dz}=⟨{\varphi }_n(z)\left|i\frac{\mathrm{\partial }}{\mathrm{\partial }z}-H(z)\right|{\varphi }_n(z)⟩.$(14)

We now turn to the particular case of a two waveguide coupler. Lai et al. [58] derive a representation of the invariant for particular Hamiltonians, including the Hamiltonian in Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) . It can be written as

(15) $H(z)=2\mathrm{\Delta }(z)K_0+\kappa (z)\left(K_{+}+K_{-}\right),$(15)

with

(16) $K_0=\frac{1}{2}\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right),K_{+}=\left(\begin{array}{cc}0& 1\\ 0& 0\end{array}\right),K_{-}=\left(\begin{array}{cc}0& 0\\ 1& 0\end{array}\right).$(16)

Since $[K_0,K_{\pm }]=\pm K_{\pm }$ and $[K_{+},K_{-}]=2K_0$, it satisfies the conditions of Lai et al. [58], and the invariant is given by

(17) $I(z)=R(z)K_{0}R(z{)}^{†},$(17)

where

(18) $R(z)=exp\left(-\frac{{\gamma }_{LR}(z)}{2}(K_{+}{e}^{i{\beta }_{LR}(z)}-K_{-}{e}^{-i{\beta }_{LR}(z)})\right).$(18)

Here ${\gamma }_{LR}(z)$ and ${\beta }_{LR}(z)$ are auxiliary real functions that parametrise the possible invariants. The following development closely follows that of Chen et al. [21] and was later applied to waveguides by Tseng [23]. Applying this procedure yields

(19) $I=\frac{1}{2}\left(\begin{array}{cc}cos{\gamma }_{LR}& e^{i{\beta }_{LR}}sin{\gamma }_{LR}\\ e^{-i{\beta }_{LR}}sin{\gamma }_{LR}& -cos{\gamma }_{LR}\end{array}\right).$(19)

Equation (19)(19) $I=\frac{1}{2}\left(\begin{array}{cc}cos{\gamma }_{LR}& e^{i{\beta }_{LR}}sin{\gamma }_{LR}\\ e^{-i{\beta }_{LR}}sin{\gamma }_{LR}& -cos{\gamma }_{LR}\end{array}\right).$(19) has Eigenvalues ${\lambda }_{\pm }=\pm 1/2$, with Eigenstates

(20) $|{\varphi }_{+}⟩=\left(\begin{array}{c}{e}^{i{\beta }_{LR}}cos\frac{{\gamma }_{LR}}{2}\\ sin\frac{{\gamma }_{LR}}{2}\end{array}\right),|{\varphi }_{-}⟩=\left(\begin{array}{c}sin\frac{{\gamma }_{LR}}{2}\\ -e^{-i{\beta }_{LR}}cos\frac{{\gamma }_{LR}}{2}\end{array}\right).$(20)

Substituting Equation (19)(19) $I=\frac{1}{2}\left(\begin{array}{cc}cos{\gamma }_{LR}& e^{i{\beta }_{LR}}sin{\gamma }_{LR}\\ e^{-i{\beta }_{LR}}sin{\gamma }_{LR}& -cos{\gamma }_{LR}\end{array}\right).$(19) into the invariance condition Equation (11)(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11) yields a coupled system of differential equations

(21) $\begin{array}{cc}\frac{d{\gamma }_{LR}}{dz}& =2\kappa sin{\beta }_{LR}\\ \frac{d{\beta }_{LR}}{dz}& =2\kappa cos{\beta }_{LR}cot{\gamma }_{LR}-2\mathrm{\Delta }.\end{array}$(21)

We note that the expression here may differ from those in other publications by a factor $2$, which arises from different definitions of $H$. The method of Lewis-Riesenfeld is indirect, in that, once ${\gamma }_{LR}$ and ${\beta }_{LR}$ are known, $\mathrm{\Delta }$ and $\kappa$ can be found. The functions ${\gamma }_{LR}$ and ${\beta }_{LR}$ are arbitrary except for their boundary conditions at $z=0$ and $z=\mathrm{\ell }$. At the ends of the coupler, there should be no coupling and thus $\kappa =0$. This gives

(22) ${\left.\frac{d{\gamma }_{LR}}{dz}\right|}_{z=0}={\left.\frac{d{\gamma }_{LR}}{dz}\right|}_{z=\mathrm{\ell }}.$(22)

This also ensures that $[I(0),H(0)]=[I(\mathrm{\ell }),H(\mathrm{\ell })]=0$, implying shared Eigenstates between $H$ and $I$. From the Eigenstates in Equation (20)(20) $|{\varphi }_{+}⟩=\left(\begin{array}{c}{e}^{i{\beta }_{LR}}cos\frac{{\gamma }_{LR}}{2}\\ sin\frac{{\gamma }_{LR}}{2}\end{array}\right),|{\varphi }_{-}⟩=\left(\begin{array}{c}sin\frac{{\gamma }_{LR}}{2}\\ -e^{-i{\beta }_{LR}}cos\frac{{\gamma }_{LR}}{2}\end{array}\right).$(20) , we observe that to start and finish in different pure states (i.e. to transfer power from one waveguide to the other) requires

(23) ${\gamma }_{LR}(0)=\pi ,{\gamma }_{LR}(\mathrm{\ell })=0,$(23)

or vice versa without loss of generality. Equations (22)(22) ${\left.\frac{d{\gamma }_{LR}}{dz}\right|}_{z=0}={\left.\frac{d{\gamma }_{LR}}{dz}\right|}_{z=\mathrm{\ell }}.$(22) and (23(23) ${\gamma }_{LR}(0)=\pi ,{\gamma }_{LR}(\mathrm{\ell })=0,$(23) ) are the boundary conditions that are necessary for a valid coupler that transfers power from one waveguide to the other. To validate the performance at all device lengths, Chen et al. [21] chose length-invariant boundary conditions of the form

(24) ${\left.\frac{d{\beta }_{LR}}{dz}\right|}_{z=0}={\left.-\frac{d{\beta }_{LR}}{dz}\right|}_{z=\mathrm{\ell }}=\frac{C}{\mathrm{\ell }},$(24)

with $C>0$ a constant. From a dimensional analysis perspective, we observe that this is the only dependence on $\mathrm{\ell }$ that leads to identical behaviour at all device lengths. Furthermore, in the interests of keeping the coupler realistic by minimising the coupling coefficient $\kappa$, Chen et al. [21] chose boundary conditions on ${\beta }_{LR}$ to maximise $sin{\beta }_{LR}$ in Equation (21)(21) $\begin{array}{cc}\frac{d{\gamma }_{LR}}{dz}& =2\kappa sin{\beta }_{LR}\\ \frac{d{\beta }_{LR}}{dz}& =2\kappa cos{\beta }_{LR}cot{\gamma }_{LR}-2\mathrm{\Delta }.\end{array}$(21)

(25) ${\beta }_{LR}(0)={\beta }_{LR}(\mathrm{\ell })=-\pi /2.$(25)

One choice that satisfies the requirements in Equations (22)(22) ${\left.\frac{d{\gamma }_{LR}}{dz}\right|}_{z=0}={\left.\frac{d{\gamma }_{LR}}{dz}\right|}_{z=\mathrm{\ell }}.$(22) –(25(25) ${\beta }_{LR}(0)={\beta }_{LR}(\mathrm{\ell })=-\pi /2.$(25) ) are third-order polynomials for ${\gamma }_{LR}(z)$ and ${\beta }_{LR}(z)$ [21].

The parameterization from $\{\kappa ,\mathrm{\Delta }\}$ to $\{{\gamma }_{LR},{\beta }_{LR}\}$ may provide additional other benefits for waveguide couplers. Tseng [23,24], leveraging the work from quantum population transfer of Ruschhaupt et al. [59], uses these parameters in the Lewis-Riesenfeld phase ${\alpha }_n$ to describe the robustness of the device to constant deviations in Hamiltonian elements. The error sensitivity as a function of invariant parameters was optimized, producing couplers that are resilient to changes in either $\kappa$ or $\mathrm{\Delta }$. This is discussed further in Section 5.

3.2. Counter-diabatic method

The method proposed by Demirplak and Rice [33] approaches the problem from a different perspective. The central observation of this approach is that an evolving non-adiabatic Hamiltonian $H_0(z)$ leads to transitions between states. These transitions arise from off-diagonal terms in the Hamiltonian, and can be countered by adding a counter-diabatic term $H_{\mathrm{C}\mathrm{D}}(z)$ to the Hamiltonian. Under the dynamics of the new, ‘corrected’ Hamiltonian $H(z)=H_0(z)+H_{\mathrm{C}\mathrm{D}}(z)$ the system remains in the initial Eigenstate of the non-corrected Hamiltonian [57], as shown in Figure 4.

Figure 4. Fractional power in the unwanted Eigenmode versus distance for a coupler designed with the counter-diabatic method. Dashed curve: power variation in the supermodes of the uncorrected Hamiltonian $H_0(z)$. Solid black curve: evolution of the corrected Hamiltonian $H_0(z)+H_{\mathrm{C}\mathrm{D}}(z)$ projected onto its Eigenstates. Red curve: same evolution, but projected onto the Eigenstates of the original Hamiltonian $\left|n(z)\right.⟩$

To find the appropriate counter-diabatic Hamiltonian, we transform the ket $\left|G(z)\right.⟩$ appearing in Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) to a basis set that evolves together with the Eigenstates of the uncorrected Hamiltonian $H_0(z)$. These Eigenstates $\left|n(z)\right.⟩$ satisfy

(26) $H_0(z)\left|n(z)\right.⟩=E_n(z)\left|n(z)\right.⟩.$(26)

with Eigenvalues $E_n(z)$. A basis transformation must be unitary, and can generally be written as

(27) $\left|\tilde{\psi }(z)\right.⟩=U_{\mathrm{C}\mathrm{D}}(z)\left|G(z)\right.⟩,$(27)

and the evolution of the transformed system is then given by the Schrödinger Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2)

(28) $\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|\tilde{\psi }(z)\right.⟩=-iH(z)\left|\tilde{\psi }(z)\right.⟩,$(28)

for a new Hamiltonian $H(z)$. To avoid transitions, this transformed Hamiltonian should be diagonal in the basis of Eigenstates $\left|n(0)\right.⟩$. $H(z)$ depends on the particular choice of $U_{\mathrm{C}\mathrm{D}}(z)$ via

(29) $H(z)=U_{\mathrm{C}\mathrm{D}}(z{)}^{†}{H}_0(z)U_{\mathrm{C}\mathrm{D}}(z)-i{U}_{\mathrm{C}\mathrm{D}}(z)\frac{\mathrm{\partial }}{\mathrm{\partial }z}{U}_{\mathrm{C}\mathrm{D}}(z{)}^{†}+U_{\mathrm{C}\mathrm{D}}(z)H_{\mathrm{C}\mathrm{D}}(z)U_{\mathrm{C}\mathrm{D}}(z{)}^{†}.$(29)

The first term can be made diagonal by choosing $U_{\mathrm{C}\mathrm{D}}(z)$ such that its row vectors are the transposed Eigenstates $\left.⟨n(z)\right|$, and so

(30) $U_{\mathrm{C}\mathrm{D}}(z)=\sum _n\left|n(0)\right.⟩\left.⟨n(z)\right|.$(30)

The remaining off-diagonal terms can be made to vanish by choosing

(31) $H_{\mathrm{C}\mathrm{D}}(z)=i\frac{\mathrm{\partial }{U}_{\mathrm{C}\mathrm{D}}(z{)}^{†}}{\mathrm{\partial }z}{U}_{\mathrm{C}D}(z)=i\sum _n\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)\left.⟨n(z)\right|.$(31)

For any given non-adiabatic Hamiltonian, we therefore have

(32) $H(z)=H_0(z)+H_{\mathrm{C}\mathrm{D}}(z)=\sum _n\left|n(z)\right.⟩E_n(z)\left.⟨n(z)\right|+i\sum _n\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)\left.⟨n(z)\right|,$(32)

where the first (diagonal) term evolves the states to match the Eigenstates at position $z$, and $H_{\mathrm{C}\mathrm{D}}(z)$ cancels out transitions resulting from the evolution. For a reference Hamiltonian $H_0(z)$ (Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) ) we explicitly find that Equation (32)(32) $H(z)=H_0(z)+H_{\mathrm{C}\mathrm{D}}(z)=\sum _n\left|n(z)\right.⟩E_n(z)\left.⟨n(z)\right|+i\sum _n\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)\left.⟨n(z)\right|,$(32) reduces to a Hamiltonian of the form

(33) $H(z)=\left(\begin{array}{cc}\mathrm{\Delta }& \kappa +\frac{i}{2}\frac{d\theta }{dz}\\ \kappa -\frac{i}{2}\frac{d\theta }{dz}& -\mathrm{\Delta }\end{array}\right),$(33)

where $\theta$ is defined as in Equation (4)(4) $tan\theta =\kappa /\mathrm{\Delta }.$(4) . In this form, we see that the counter-diabatic Hamiltonian intuitively manifests as a correction to the coupling term $\kappa$, and this correction is related to the adiabaticity criterion $\eta$ (Equation (7)(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) ).

3.3. Transitionless-tracking algorithm

Berry’s transitionless tracking algorithm [34] takes a similar approach to suppressing transitions between states. Instead of changing the basis of the Hamiltonian to one in which the elements are diagonal in the basis of initial states, the approach seeks to find a unitary operator that preserves adiabatic states under evolution. These adiabatic states correspond to the states traversed by the system if the change in the Hamiltonian is performed adiabatically, and are, to within a phase factor, the same as the instantaneous Eigenstates $\left|n(z)\right.⟩$

(34) $\left|{\psi }_n(z)〉=e^{i{\varphi }_n(z)}\right|n(z)\text{ rightrangle}$(34)

where ${\varphi }_n(z)$ incorporates the integrated phase accumulation from the energy $E_n(z)$ of the state, together with the geometric phase. A unitary transformation that preserves the probability distributions of states under evolution from $0$ to $z$ is

(35) $U_{\mathrm{T}\mathrm{T}}(z)=\sum _n{e}^{i{\varphi }_n(z)}\left|n(z)\right.⟩\left.⟨n(0)\right|.$(35)

The Hamiltonian that performs such a transformation can then be found. Substituting the relation $\left|G(z)\right.⟩=U_{\mathrm{T}\mathrm{T}}(z)\left|\psi (0)\right.⟩$ into Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) we find

(36) $\begin{array}{rl}& H(z)=i\frac{\mathrm{\partial }{U}_{\mathrm{T}\mathrm{T}}(z)}{\mathrm{\partial }z}{U}_{\mathrm{T}\mathrm{T}}(z{)}^{†}=\sum _n\left|n(z)\right.⟩E_n(z)⟨n(z)|\\ & +i\sum _n\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)⟨n(z)|-i\sum _n\left|n(z)\right.⟩⟨n(z)|\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)⟨n(z)|\\ & \equiv H_0(z)+H_{\mathrm{T}\mathrm{T}}(z),\end{array}$(36)

where $H_0(z)$ is the original Hamiltonian and $H_{\mathrm{T}\mathrm{T}}(z)$ is a correction that prevents transitions. In the particular case of a two waveguide coupler with $H_0(z)$ as in (2), the transitionless tracking algorithm yields the same correction as the counter-diabatic approach (33).

3.4. Unitary transformation enabling photonics applications

Both the counter-diabatic and transitionless tracking methods provide a correction term to the Hamiltonian. These corrections are generally complex-valued, and indeed the corrections are purely imaginary in the case of a two waveguide coupler. In the context of atomic physics and quantum state driving, their implementation is straightforward, and the imaginary component is the phase of the driving field. In optics however, the total Hamiltonian must be real in order to correspond to a passive optical device. Following Tseng [48], this may be achieved through one final unitary transformation $V_{\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{e}}=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(exp(i\chi /2),exp(-i\chi /2))$ which applies a rotation in the real-imaginary plane. The rotation angle $\chi$ is real and yet to be determined. Applying this unitary transformation to the Schrödinger equation, as before (see Equation (28)(28) $\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|\tilde{\psi }(z)\right.⟩=-iH(z)\left|\tilde{\psi }(z)\right.⟩,$(28) ), we find the transformed Hamiltonian

(37) $H_{V\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}}(z)=V_{\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{e}}(z)H(z)V_{\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{e}}^{†}(z)-i{V}_{\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{e}}(z)\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}_{\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{e}}^{†}(z).$(37)

If we apply this to the explicit form of the counter-diabatic Hamiltonian Equation (33)(33) $H(z)=\left(\begin{array}{cc}\mathrm{\Delta }& \kappa +\frac{i}{2}\frac{d\theta }{dz}\\ \kappa -\frac{i}{2}\frac{d\theta }{dz}& -\mathrm{\Delta }\end{array}\right),$(33) , we find that off-diagonal elements take the form $(\kappa \mp \frac{i}{2}(d\theta /dz))exp(\pm i\chi )$. By choosing $\mathrm{t}\mathrm{a}\mathrm{n}\chi =\frac{1}{2\kappa }\frac{d\theta }{dz}$ these elements are real, and the transformed Hamiltonian is then

(38) $H_{V\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}}(z)=\left(\begin{array}{cc}\mathrm{\Delta }-\frac{1}{2}\frac{d}{dz}arctan\left(\frac{1}{2\kappa }\frac{d\theta }{dz}\right)& \sqrt{{\kappa }^2+{(\frac{1}{2}\frac{d\theta }{dz})}^2}\\ \sqrt{{\kappa }^2+{(\frac{1}{2}\frac{d\theta }{dz})}^2}& -\mathrm{\Delta }+\frac{1}{2}\frac{d}{dz}arctan\left(\frac{1}{2\kappa }\frac{d\theta }{dz}\right)\end{array}\right).$(38)

Thus, for this $\chi$, $H_{V\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}}(z)$ corresponds to a real matrix which can be implemented as a passive waveguide coupler. The same unitary transformation can be used in the transitionless tracking approach. Indeed, the counter-diabatic and transitionless tracking approach are closely related.

3.5. Equivalence between transitionless tracking and the counter-diabatic approach

There is a direct connection between the transitionless tracking algorithm (Section 3.3) and the counter-diabatic approaches (Section 3.2). As an overall comment, the transitionless tracking algorithm preserves adiabatic states, whereas the counter-diabatic approach preserves the Eigenstates $\left|n(z)\right.⟩$ of the instantaneous Hamiltonian. Because these are the same except for a phase factor, the two correcting Hamiltonians must be the same except along the diagonal. This can be seen directly by examining the diagonal components of $H_{\mathrm{T}\mathrm{T}}(z)$

(39) $\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}\left[H_{\mathrm{T}\mathrm{T}}(z)\right]=i\sum _n\left|n(z)\right.⟩\left.⟨n(z)\right|\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)\left.⟨n(z)\right|,$(39)

and so, comparing with Equation(36), we find

(40) $H_{\mathrm{C}\mathrm{D}}(z)=H_{\mathrm{T}\mathrm{T}}(z)-\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}\left[H_{\mathrm{T}\mathrm{T}}(z)\right].$(40)

The corrected Hamiltonian $H_{\mathrm{C}\mathrm{D}}(z)$ from the counter-diabatic correction is therefore equal to the purely off-diagonal part of the transitionless tracking algorithm $H_{\mathrm{T}\mathrm{T}}(z)$. The two corrections therefore suppress the transitions between Eigenstates, but these Eigenstates may still undergo changes in phase. The transitionless tracking algorithm ensures that these phases evolve in the same way as those of pure adiabatic states.

Berry [34] notes that the phases in Equation (34)(34) $\left|{\psi }_n(z)〉=e^{i{\varphi }_n(z)}\right|n(z)\text{ rightrangle}$(34) can be chosen freely; doing so preserves the probability distribution amongst the states, although the states themselves evolve with a different phase than that expected from adiabatic evolution. Choosing $\varphi (z)=0$ leads directly to the counter-diabatic approach.

3.6. Equivalence of transitionless tracking and the invariant approach

Similarly, we find that a Hamiltonian computed by the transitionless tracking method can also be generated through an invariant approach. First, consider a coupler designed by transitionless tracking. The Hamiltonian always drives states adiabatically (illustrated for a particular case in Section 4). Chen et al. [21] propose that an operator of the form

(41) $I(z)=\sum _n\left|n(z)\right.⟩{\lambda }_n\left.⟨n(z)\right|,$(41)

satisfies the invariance condition Equation (11)(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11) with different choices of ${\lambda }_n$. This invariant makes intuitive sense in the context of Figure 4, since the basis of the original Hamiltonian $\left|n(z)\right.⟩$ remains stationary under the dynamics of $H_0(z)+H_{\mathrm{C}D}(z)$. The simple choice ${\lambda }_n=E_n(0)$ ensures that the invariant and original Hamiltonian are identical at $z=0$, giving the invariant a more grounded interpretation. In Appendix A1 we show explicitly that $I(z)$ satisfies Equation (11)(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11) . The result also demonstrates that the choice of Eigenvalues of the invariant ${\lambda }_n$ hold no great significance in terms of the final coupler produced. In the opposite direction, we can ask if, for every Hermitian invariant $I(z)$ satisfying Equation (11)(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11) , there exists a Hamiltonian $H_0(z)$, which, when applying the transitionless tracking method, gives the same $H(z)$ as the invariant method. Chen et al. [21] propose that the choice

(42) $H_0(z)=\sum _n\left|n(z)\right.⟩\left(i\left.⟨n(z)\right|\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩-\frac{d{\alpha }_n}{dz}\right)\left.⟨n(z)\right|,$(42)

accomplishes this, where $\left|n(z)\right.⟩$ again represents the Eigenvalues of the invariant. The Lewis-Riesenfeld phase ${\alpha }_n(z)$ can be calculated from the invariant approach Equation (14)(14) $\frac{d{\alpha }_n}{dz}=⟨{\varphi }_n(z)\left|i\frac{\mathrm{\partial }}{\mathrm{\partial }z}-H(z)\right|{\varphi }_n(z)⟩.$(14) , and it replaces the Berry phase ${\phi }_n(z)$.

Unitary transforms, however, result in a need to modify this: if $H(z)$ is real when produced by the transitionless tracking algorithm, and $I(z)$ from Equation (41)(41) $I(z)=\sum _n\left|n(z)\right.⟩{\lambda }_n\left.⟨n(z)\right|,$(41) is real since $\left|n(z)\right.⟩$ are real (from a real $H_0(z)$), then Equation (11)(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11) is not satisfied since one term is real and the other imaginary. Hence, we propose a different invariant

(43) $I_V=\sum _{n}V(z)\left|n(z)\right.⟩{\lambda }_n⟨n(z)|V^{†}(z),$(43)

with $V$ as some unitary transform and $\left|n(z)\right.⟩$ is the basis state for the original Hamiltonian $H_0(z)$. In Appendix A2 we explicitly verify that this satisfies Equation (11) with $H(z)=H_v(z)$ from Equation (37)(37) $H_{V\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}}(z)=V_{\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{e}}(z)H(z)V_{\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{e}}^{†}(z)-i{V}_{\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{e}}(z)\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}_{\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{e}}^{†}(z).$(37) . The result is applicable to any unitary transformation $V$, including $V_{\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{e}}$ in Section 3.4. It is possible to generate a plot similar to Figure 4 for the unitary transformed dynamics, where projecting the state vector $|G⟩$ onto the transformed basis $V\left|n(z)\right.⟩$ shows zero changes in mode power. Again, this makes intuitive sense—stationary states should become the basis of the invariant. In the context of waveguide couplers, however, this basis not only no longer exists but also it is unphysical since $V\left|n(z)\right.⟩\in {\mathbb{C}}^2$. In the other direction of equivalence, one can obtain an initial Hamiltonian $H_0(z)$ by transforming $\left|n(z)\right.⟩$ to $V\left|n(z)\right.⟩$ in Equation (42)(42) $H_0(z)=\sum _n\left|n(z)\right.⟩\left(i\left.⟨n(z)\right|\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩-\frac{d{\alpha }_n}{dz}\right)\left.⟨n(z)\right|,$(42) .

Thus, all three of the presented methods are equivalent, at the operator level, even when constraining the Hamiltonian through unitary transformation to produce physical waveguides.

4. Adiabatic processes in realistic waveguides

So far, we have reviewed adiabatic processes on the basis of coupled mode theory (CMT) using Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) as a starting point, with the parameters $\kappa (z)$ and $\mathrm{\Delta }(z)$ characterizing the coupler and its performance (e.g. the crosstalk $\mu$). While this approach provides a wealth of information on coupler properties, the conceptual link to physical devices is less obvious. In this section, we review pathways and implications of implementing adiabatic couplers in practice. The modes and full propagation properties of arbitrary geometries can be calculated with commercial numerical mode solvers (e.g. COMSOL, Lumerical, etc.), and can be used as a tool for evaluating the coupled mode formalism. In this Review, we show results generated by finite element method (FEM) calculations implemented in COMSOL.

The link between coupled mode theory and physical waveguides can be approached in two ways: (i) a physical embodiment of a coupler can be used to calculate $\kappa (z)$ and $\mathrm{\Delta }(z)$ to solve Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) and predict coupler properties; (ii) vice versa, a desired distribution of $\kappa (z)$ and $\mathrm{\Delta }(z)$ targeting certain device characteristics can be given a physical embodiment. We now review both cases, comparing the crosstalk predicted by various methods with full finite element calculations. Section 4.1 reviews the properties of a previously reported two-waveguide adiabatic coupler [47] as an illustrative example and as starting point for the subsequent discussion. Section 4.2 then presents typical approaches for applying two different protocols to this example physical system: the counter-diabatic correction of Sec. 3.4 in Section 4.2.1, and the FAQUAD approach from Sec. 2 in Section 4.2.2. Section 4.2.3 concludes with a comparison of the physical embodiment of the $n=\pm 2$ polynomials of Figure 3 with a reasonable choice of $\mathrm{\Gamma }$.

Below we discuss the properties of example couplers and their scaled counterparts. The scaled device is related to the original device by stretching or compressing without other modifications. As we showed in Section 3.2 the CD correction guarantees vanishing crosstalk but requires a different correction for every choice of length. In the spirit of analyzing a device’s crosstalk envelope however, below we apply the CD correction once for a single length, and present how the crosstalk varies as the resulting device is scaled.

Figure 5. Example adiabatic device from Ng et al. [47]. (a) Finite element calculation and definition of physical parameters $w_1$, $w_2$, and $s$. (b) Resulting effective index of the isolated and coupled waveguides. (c) Computed $\mathrm{\Delta }$, $\mathrm{\Gamma }$ and $\kappa$, and (d) $\theta$. (e) Resulting $\eta$ (dashed line: $\eta =1$.) (f) Crosstalk according to Louisell (Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) ), coupled mode theory (CMT, solution of Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) ), and finite element calculations

4.1. Parabolic coupler

As discussed in Section 2, the coupled mode formalism is agnostic to the waveguides’ physical embodiment: two parameters $\kappa$ and $\mathrm{\Delta }$ describe infinitely many materials, wavelengths and geometries, which in a chip-based, two-waveguide scenario can include the width and height of each waveguide, and their relative separations. In order to simplify the discussion and to reduce the physical parameter space, we consider a previously reported example [47], consisting of two 1D waveguides (WGs) at $\lambda =1.55\mu \mathrm{m}$, formed by silica ($n_{\mathrm{S}\mathrm{i}{\mathrm{O}}_2}=1.444$ [60]) in air, under transverse magnetic (TM) polarization, as shown in Figure 5(a). In this case, the only relevant physical parameters are the widths of WG1 and WG2 ($w_1(z)$ and $w_2(z)$, respectively), and their edge-to-edge spacing $s(z)$. This representative device is a judicious embodiment of the Figure 1(b) schematic; the width of the input WG1 decreases linearly from 700 nm to 300 nm, and the width of WG2 increases from 300 nm to 700 nm, so that $\mathrm{\Delta }$ is approximately a straight line with negative slope with $\mathrm{\Delta }=0$ at the centre of the device; the separation varies parabolically such that $s=4\mu \mathrm{m}$ at the edges $s=500\mathrm{n}\mathrm{m}$ in the centre, ensuring that the coupling coefficient $\kappa$ vanishes at the edges and peaks in the centre. Combining these characteristics broadly satisfies typical requirements for $\mathrm{\Delta }$ and $\kappa$ of an adiabatic coupler [46]. The white curves in Figure 5(a) give a geometrical outline of the coupler. The colors show the Poynting vector for the particular example of $\mathrm{\ell }=40\mu \mathrm{m}$, and confirm that most of the light couples between the waveguides in the region where they approach each other closely.

Figure 5(b) shows the calculated effective index $n_{\mathrm{e}\mathrm{f}\mathrm{f},\mathrm{j}},(j=1,2)$ of the Eigenmodes of the full two-waveguide system versus position. The dashed lines in Figure 5(b) show the corresponding effective index ${\tilde{n}}_{\mathrm{e}\mathrm{f}\mathrm{f},\mathrm{j}},(j=1,2)$ of the modes of the isolated waveguides. The parameters in Figure 5(b) yield the parameters $\kappa$, $\mathrm{\Delta }$ and $\mathrm{\Gamma }$, shown in Figure 5(c) which are obtained directly from

(44) $\begin{array}{rl}& \mathrm{\Gamma }=k_0(n_{\mathrm{e}\mathrm{f}\mathrm{f},1}-n_{\mathrm{e}\mathrm{f}\mathrm{f},2})/2,\\ & \mathrm{\Delta }=k_0({\tilde{n}}_{\mathrm{e}\mathrm{f}\mathrm{f},1}-{\tilde{n}}_{\mathrm{e}\mathrm{f}\mathrm{f},2})/2,\\ & \kappa =\sqrt{{\mathrm{\Gamma }}^2-{\mathrm{\Delta }}^2}.\end{array}$(44)

These can be used to obtain $\theta$ from Equation (4)(4) $tan\theta =\kappa /\mathrm{\Delta }.$(4) (see Figure 5(d)). Subsequently, $\eta$ can be obtained from Equation (7)(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) , shown in Figure 5(e) for different device lengths $\mathrm{\ell }$. Recall that the adiabatic condition is $\eta \ll 1$: this condition is satisfied at all points along the device for $\mathrm{\ell }\stackrel{>}{\approx \phantom{{ }_2}}50\mu \mathrm{m}$ (green curve in Figure 5(e)), and is satisfied near the edges when $\mathrm{\ell }=15\mu \mathrm{m}$ (blue curve in Figure 5(e)). Finally, $\mathrm{\ell }\sim 27.5\mu \mathrm{m}$ represents an intermediate case where $\eta \le 1$ along the entire device.

Knowing both physical parameters and coupled-mode parameters, the crosstalk $\mu$ can be calculated in three different ways. The results are summarized in Figure 5(f), which shows $\mu$ versus device length $\mathrm{\ell }$. The blue curve gives the result of integrating Louisell’s result (Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) ) [46], and the orange curve gives the result obtained by numerically integrating the coupled mode equations (Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) ). The agreement between these curves shows that Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) is accurate for the parameters in this example, whereby the largest discrepancy occurs at short lengths where Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) is not expected to be valid [46]. The black curve shows $\mu$ obtained from a full two-dimensional finite element calculation of the physical system in Figure 5(a), showing excellent agreement with the previous methods, confirming the validity of all methods to describe this physical system.

4.2. Application of shortcuts and protocols

With knowledge of all important parameters for our chosen system, we are now in a position to apply the procedures presented in Section 3, and review their impact on both performance and physical embodiment of the device. Note that there are still many combinations of $w_1$, $w_2$, and $s$ that lead to the same $\kappa$ and $\mathrm{\Delta }$, even keeping the material and wavelength constant as we do here. Therefore, when choosing the physical representation for a given set of parameters, we fix the bottom waveguide in its current configuration (i.e. $w_1$ decreases linearly from 700 nm to 300 nm), and modify the edge-to-edge separation $s$ and top waveguide width $w_2$ to yield a desired $\mathrm{\Delta }$ and $\kappa$ at each $z$. For a fixed $w_1$, there is a unique value of $w_2$ that provides a desired $\mathrm{\Delta }$, which describes the uncoupled waveguides. Once this is known, there is a unique value of $s$ that yields the desired $\kappa$. This approach is valid for waveguides that are not too strongly coupled, such that the off-diagonal- terms in Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) are equal [61,62]. Our choice significantly reduces the potential parameter space, but is representative of typical experimental approaches [63], and enables a direct comparison between the physical embodiment and coupled mode theory, providing insight into the applicability of this formalism for practical applications.

4.2.1. Counter-diabatic correction

Since the counter-diabatic, transitionless tracking and invariant approaches are all equivalent, we may consider the performance of a counter-diabatic device to be indicative of all three. The counter-diabatic correction for a device of length $\mathrm{\ell }$ is obtained by replacing $H$ in Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) with $H_V$ in Equation (38)(38) $H_{V\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}}(z)=\left(\begin{array}{cc}\mathrm{\Delta }-\frac{1}{2}\frac{d}{dz}arctan\left(\frac{1}{2\kappa }\frac{d\theta }{dz}\right)& \sqrt{{\kappa }^2+{(\frac{1}{2}\frac{d\theta }{dz})}^2}\\ \sqrt{{\kappa }^2+{(\frac{1}{2}\frac{d\theta }{dz})}^2}& -\mathrm{\Delta }+\frac{1}{2}\frac{d}{dz}arctan\left(\frac{1}{2\kappa }\frac{d\theta }{dz}\right)\end{array}\right).$(38) , as outlined in Ref [48]. Since $z=\zeta \mathrm{\ell }$, it follows that $H_V(z/\mathrm{\ell })$ is different for every choice of $\mathrm{\ell }$. As an illustrative example, we apply the CD correction to the waveguide system of Figure 5 via the parameters in Figure 5(c) for a device of length $\mathrm{\ell }=20\mu \mathrm{m}$. The results are shown in Figure 6(a): relative to the original structure, the CD protocol induces a strong modulation of the corrected $\mathrm{\Delta }$ and $\kappa$ in the central region, which overall increases $\mathrm{\Gamma }$. The parameters of the physical embodiment, given our constraints on the bottom waveguide, are shown in Figure 6(b): $w_2$ and $s$ vary in the centre, leading to an increase in the coupling constant $\kappa$ in the center of the device. The strongly oscillating $\theta$ (blue curve in Figure 6(c)), leads to $\eta >1$ for our value of $\mathrm{\ell }$, indicating that this CD-corrected device is less adiabatic than the original. Calculating the crosstalk $\mu$ using different methods, shown in Figure 6(d), reveals that while this design gives the desired low crosstalk $\mu$ for this specific $\mathrm{\ell }$ (black arrow in Figure 6(d)), the scaling of this profile for increasing lengths leads to an overall larger envelope than the original structure. Finally, while the finite element calculations of the physical device faithfully reproduce the overall crosstalk envelope, shown in Figure 6(d), as well as the salient features predicted by Equations (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) and (6(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) ), the physical device has a larger crosstalk at the design length $\mathrm{\ell }$ than Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) predicts, as a result of the large curvature of the top waveguide over wavelength-scale distances, leading to the leaking of radiation into the bottom waveguide.

Figure 6. CD (left) and FAQUAD (right) corrections to the device in Figure 5. (a) $\mathrm{\Delta }$, $\mathrm{\Gamma }$ and $\kappa$ obtained by applying the CD correction (Equation (38)(38) $H_{V\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{l}}(z)=\left(\begin{array}{cc}\mathrm{\Delta }-\frac{1}{2}\frac{d}{dz}arctan\left(\frac{1}{2\kappa }\frac{d\theta }{dz}\right)& \sqrt{{\kappa }^2+{(\frac{1}{2}\frac{d\theta }{dz})}^2}\\ \sqrt{{\kappa }^2+{(\frac{1}{2}\frac{d\theta }{dz})}^2}& -\mathrm{\Delta }+\frac{1}{2}\frac{d}{dz}arctan\left(\frac{1}{2\kappa }\frac{d\theta }{dz}\right)\end{array}\right).$(38) ) to the parameters in Figure 5(c),(d) for $\mathrm{\ell }=20\mu m$, and (b) resulting physical parameters defined in Figure 5(a). (c) Comparison of $\theta /\pi$ (left axis) and $\eta$. (d) Resulting crosstalk according to Louisell (Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) , blue), coupled mode theory (solution of Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) , orange), and FEM calculations (black). (e) obtained by applying the FAQUAD protocol to the parameters in Figure 5(c), and (d), and (f), resulting physical parameters as defined in Fig. 5(a). (g) Comparison of $\theta /\pi$ and $\eta$ at $\mathrm{\ell }=10\mu m$: $\eta$ is constant as required. (h) Resulting crosstalk according to Louisell (Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) , dashed blue), coupled mode theory (CMT, solution of Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) , orange), and FEM calculations (black)

4.2.2. FAQUAD correction

The fast quasi-adiabatic (FAQUAD) correction [64] is a design protocol which aims at homogeneously distributing adiabaticity during propagation, which means that $\eta$ in Equation (7)(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) is constant. While there are many ways to achieve this, a typical approach [64] integrates Equation (7)(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) using $\theta$ and $\mathrm{\Gamma }$ of a chosen device. This leads to a transformation of $z$: intuitively, Equation (7)(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) informs the degree by which $z$ should be stretched or compressed along the device so as to distribute adiabaticity equally, thereby avoiding scenarios wherein the adiabatic condition is satisfied only for certain regions. This is especially important for short devices—see, for example, the blue curve in Figure 5(e), which shows that for a short device length Equation (7)(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) is satisfied at the edges but not in the centre.

The parameters obtained by the FAQUAD correction to the device in Figure 5 are shown in Figure 6(e). The corresponding physical device parameters, shown in Figure 6(f), are consistent with the expectations from the adiabaticity criterion at short lengths: at the edges, the device already satisfies the adiabaticity criterion so that $w_2$ and $s$ vary faster, whereas in the centre these parameters vary more slowly. The resulting values of $\theta$ and $\eta$ are shown in Figure 6(f): $\eta$ is flat and constant, as required by the FAQUAD condition, whereas $\theta$ approaches a straight line. The calculated crosstalk for this device using different methods is shown in Figure 6(h). Note that the crosstalk envelope scales as $\mu \propto {\mathrm{\ell }}^{-2}$, consistent with earlier results for $n=0$ [47,56]. However, this scales less rapidly than for the original coupler, which follows $\mu \propto {\mathrm{\ell }}^{-6}$ (Figure 5(f)). This exemplifies the scenario whereby the FAQUAD procedure pushes the first zero to a shorter device length, at the expense of increasing the crosstalk envelope.

Figure 7. Polynomial approach with $\mathrm{\Gamma }=\pi /35\mu {\mathrm{m}}^{-1}$. (a) $\mathrm{\Delta }$, $\mathrm{\Gamma }$ and $\kappa$ obtained from Equation (9)(9) $P{{ }^{\mathrm{\prime }}}_n(\zeta )=\frac{(2n+1)!}{{(n!)}^2}(\zeta (1-\zeta ){)}^n,$(9) for $n=2$. (b) Retrieved physical parameters defined in Figure 5(a). (c) Comparison of $\theta /\pi$ and $\eta$ at $\mathrm{\ell }=35\mu \mathrm{m}$. (d) Resulting crosstalk according to Louisell (Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) , blue), coupled mode theory (solution of Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) , orange), and FEM calculations (black). (e) $\mathrm{\Delta }$, $\mathrm{\Gamma }$ and $\kappa$ obtained from Equation (9)(9) $P{{ }^{\mathrm{\prime }}}_n(\zeta )=\frac{(2n+1)!}{{(n!)}^2}(\zeta (1-\zeta ){)}^n,$(9) for $n=-2$. (f) Retrieved physical parameters defined in Figure 5(a). (g) Comparison of $\theta /\pi$ and $\eta$ at $\mathrm{\ell }=35\mu \mathrm{m}$. (h) Resulting crosstalk according to Louisell (Equation (6)(6) $\mu =\frac{1}{4}{\left|{\int }_0^{\mathrm{\ell }}\frac{d\theta }{dz}{e}^{-2i{\int }_0^z\mathrm{\Gamma }(z^{\mathrm{\prime }})d{z}^{\mathrm{\prime }}}dz\right|}^2=\frac{1}{4}{\left|{\int }_0^{\rho (\mathrm{\ell })}\frac{d\theta }{d\rho }{e}^{-2i\rho }d\rho \right|}^2,$(6) , blue), coupled mode theory (CMT, solution of Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) , orange), and FEM calculations (black)

4.2.3. Polynomial profiles

An alternative approach for adiabatic coupler design relies on appropriately defining the $\kappa$ and $\mathrm{\Delta }$ profiles at constant $\mathrm{\Gamma }$, via polynomials of increasing order $n$ (Equation (9)(9) $P{{ }^{\mathrm{\prime }}}_n(\zeta )=\frac{(2n+1)!}{{(n!)}^2}(\zeta (1-\zeta ){)}^n,$(9) ). We have already discussed how increasing $n$ leads to a reduced device crosstalk envelope as $\mu \propto {\mathrm{\ell }}^{-2(n+1)}$ for $n\ge 0$ and when $\mathrm{\ell }\to \mathrm{\infty }$ [47]. We now give $\theta$ and $\mathrm{\Delta }$ a physical embodiment for $n=\pm 2$. Note that, in contrast to counter-diabatic/transitionless tracking and FAQUAD, this protocol relies on designing a device on the basis of pre-determined coupled-mode parameters, rather than correcting an initial profile. In order to implement the desired structure, the only new choice is that of the parameter $\mathrm{\Gamma }=\pi /L_b$, where $L_b$ is a characteristic beat length [46]. We choose $L_b=35\mu \mathrm{m}$ as shown in the blue curve of Figure 7(a), which is comparable to the minimum $\mathrm{\Gamma }$ of our example device (blue curve in Figure 5(c)). The resulting values of $\kappa$ and $\mathrm{\Delta }$ when $n=+2$ are shown in the orange and red curves of Figure 7(a), respectively. With a linearly decreasing bottom waveguide width (dashed orange line in Figure 7(b)), the edge-to-edge separation follows again a quasi-parabolic shape (blue line in Figure 7(b)); in contrast to all cases considered so far however, the width of top waveguide $w_2$ also decreases linearly near the edges (orange line in Figure 7(b)), with a cross-over point in the centre. As earlier, this corresponds to increasing and then decreasing the coupling between waveguides (as per $\kappa$) by first bringing them together and further apart, while changing their relative widths gradually near the edges and rapidly in the centre near the cross-over point (following $\mathrm{\Delta }$). The resulting $\theta$, shown in Figure 7 (blue line), follows the expected functional form (see also Figure 3). The corresponding value of $\eta$ calculated for $\mathrm{\ell }=L_b$ is shown in Figure 7 (orange line) and suggests that, strictly speaking, the adiabatic condition of Equation (7)(7) $\frac{1}{2\mathrm{\Gamma }}\frac{d\theta }{dz}\equiv \eta \ll 1.$(7) is only weakly satisfied since $\eta \stackrel{<}{\approx \phantom{{ }_2}}1$. However, as shown in Figure 3(c), this device exhibits the most rapid drop in crosstalk ($\mu \propto {\mathrm{\ell }}^{-6}$), despite having the smallest average $\mathrm{\Gamma }$ across its length.

Finally, we consider the equivalent case where $n=-2$, exemplified by the parameters in Figure 7(e). Relative to $n=+2$, this device approaches the structures of a directional-coupler: the inter-waveguide separation varies slowly in the centre, with rapid changes only at the edges (blue curve in Figure 7(f)); the top waveguide width $w_2$ varies weakly with length (solid orange curve in Figure 7(f)), and although $\eta \ll 1$ for much of the device (Figure 7(g)), the envelope decreases slowly, with periodic variations that are reminiscent of a directional coupler – compare with Figure 1(b).

Figure 8. Bandwidth calculations. (a) CMT (solid lines) and FEM calculations (circles) for the crosstalk as a function of wavelength for the physical device shown in Figure 5(a), at $\mathrm{\ell }=20\mu \mathrm{m}$ (blue), $\mathrm{\ell }=105\mu \mathrm{m}$ (red) and $\mathrm{\ell }=150\mu \mathrm{m}$ (yellow). (b) Equivalent calculations for the device with the physical profile shown in Figure 6(b), where the CD protocol has been applied at $\lambda =1550\mathrm{n}\mathrm{m}$ and $\mathrm{\ell }=20\mu \mathrm{m}$. Dashed lines indicate the crosstalk envelope

4.2.4. Bandwidth

One key advantage of adiabatic couplers with respect to directional couplers is their larger bandwidth which, as illustrated in Figure 1, is due to slow changes in a single mode’s spatial profile, rather than interference effects between multiple modes. In practice however, the modes supported by an adiabatic device weakly couple as its spatial profile changes, so that the drop in the crosstalk envelope for increasing $\mathrm{\ell }$ is accompanied by periodic local minima due to interference effects (see, for example, Figures 1(b), 5(f) and 6(d)).

Similar behaviour occurs in the wavelength dependence. Figure 8(a) shows the crosstalk versus wavelength for the adiabatic device of Figure 5(a) for three different representative lengths, and ignoring material dispersion. We consider $\mathrm{\ell }=20\mu \mathrm{m},105\mu \mathrm{m},150\mu \mathrm{m}$, which, respectively, correspond to a short device, and a local minimum/maximum of the crosstalk envelope at $\lambda =1550\mathrm{n}\mathrm{m}$ (see also Figure 5(f)). Increasing $\mathrm{\ell }$ leads to an overall decrease in the bandwidth envelope (dashed line), but with periodic local minima due to interference. Note in particular the good agreement between FEM calculations (dots) and the CMT approach (solid lines), obtained by numerically solving Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) after separately calculating $\kappa (z)$, $\mathrm{\Delta }(z)$ and $\mathrm{\Gamma }(z)$ at every wavelength via Equation (44)(44) $\begin{array}{rl}& \mathrm{\Gamma }=k_0(n_{\mathrm{e}\mathrm{f}\mathrm{f},1}-n_{\mathrm{e}\mathrm{f}\mathrm{f},2})/2,\\ & \mathrm{\Delta }=k_0({\tilde{n}}_{\mathrm{e}\mathrm{f}\mathrm{f},1}-{\tilde{n}}_{\mathrm{e}\mathrm{f}\mathrm{f},2})/2,\\ & \kappa =\sqrt{{\mathrm{\Gamma }}^2-{\mathrm{\Delta }}^2}.\end{array}$(44) .

We next analyze the bandwidth for the device with the physical profile shown in Figure 6(b), which was obtained by applying the CD protocol to the device of Figure 5(a) at $\mathrm{\ell }=20\mu \mathrm{m}$ and $\lambda =1550\mathrm{n}\mathrm{m}$. In this case, CMT calculations show perfect transfer where the device protocol was applied, with an increase in $\mu$ at other wavelengths. Nevertheless, the crosstalk envelope at $\mathrm{\ell }=20\mu \mathrm{m}$ is reduced with respect to the original device by approximately 6 dB. This conclusion is consistent with recent experiments [41,55,65] reporting that STA devices can achieve good bandwidth at relatively shorter lengths. Physically scaling this device to longer lengths, however, shows the same overall increase in crosstalk envelope observed when comparing Figures 5(f) and 6(d). Note the good agreement between FEM and CMT calculations, except for the shortest device where the CD protocol was applied, as a result of radiation leakage into the bottom waveguide due to the sharp curvature of the top waveguide at short lengths. Although this might be improved by other choices of $w_1$, $w_2$ and $s$ which produce the same $\kappa$, $\mathrm{\Delta }$ and $\mathrm{\Gamma }$, this example suggests that there are some difficulties in translating the extreme performance of CD devices at extremely short lengths into physically realistic devices.

Finally, the analysis of the dependence of crosstalk on $\mathrm{\ell }$ and $\lambda$ highlights the importance of considering the crosstalk envelope, and more generally the envelope of any parameter of interest, since it provides a broader view of the robustness of the device to interference effects.

5. Noise sensitivity

Adiabatic couplers are expected to be less sensitive to fabrication errors than directional couplers [11,66,67]. Estimating the performance of a coupler in the presence of fabrication errors is difficult as detailed knowledge of the fabrication errors is required. Furthermore, the connection between these physical errors and coupled mode theory is non-trivial (see Section 4). However, the frequency of a noise source is conserved between the two descriptions. We therefore evaluate the device performance by simulating noise on $\kappa$, $\mathrm{\Delta }$ and determining the crosstalk by numerically solving Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) . We choose devices all with the same length $\mathrm{\ell }=1$ that exhibit perfect coupling in the absence of noise. We then apply Gaussian noise characterized by a variance ${\sigma }^2$, such that $\sigma$ is proportional to the average of $\mathrm{\Gamma }$ in each device consistent with $\mathrm{\ell }=1$. This approach, which is agnostic to a particular type of fabrication error, allows us to compare the properties of different types of designs discussed in Section 4 in including one that was optimized for insensitivity to noise [23].

We take the noise to be Gaussian, characterized by a variance $\sigma$ and with an exponentially decaying correlation with characteristic length ${\mathrm{\ell }}_c$, which can be much smaller, similar or much longer than the coupler length $\mathrm{\ell }$. Our implementation follows that of Farahmand and de Sterke [68], according to whom the deviation $\delta \mathrm{\Psi }$ (where $\mathrm{\Psi }$ represents either $\kappa$ or $\mathrm{\Delta }$) at $z=z_1$ can be found from that at $z=z_0$, via the Gaussian distribution

(45) $f(\delta \mathrm{\Psi }(z_1)|\delta \mathrm{\Psi }(z_0))=\frac{1}{\sigma \sqrt{2\pi (1-r^2)}}exp\left(-\frac{1}{2{\sigma }^2(1-r^2)}{(\delta \mathrm{\Psi }(z_1)-r\delta \mathrm{\Psi }(z_0))}^2\right),$(45)

where $r(z)=exp(-|z|/{\mathrm{\ell }}_c)$. Steps are taken in $z$ in intervals of $0.004{\mathrm{\ell }}_c$ or $0.001\mathrm{\ell }$, whichever is smaller. This model is Gaussian to all orders, with $\sqrt{1-r^2}$ scaling the standard deviation $\sigma$ for each step. Since shorter couplers are in general more sensitive and have higher $\mathrm{\Gamma }$, we choose a relative error $\epsilon =1\mathrm{\%}$ with respect to the average value of $\mathrm{\Gamma }$

(46) $\sigma =\frac{\epsilon }{\mathrm{\ell }}{\int }_0^{\mathrm{\ell }}\mathrm{\Gamma }(z)dz.$(46)

Figure 9 (a) and (f) illustrate the resulting variation of parameters $\kappa$ and $\mathrm{\Delta }$ from the ideal along the length of the device for two different correlation lengths. Note that $\epsilon$ has been exaggerated ($\epsilon$=10%) to clearly show how the correlation length affects the noise.

Figure 9. Top row solid curves: exaggerated noise model with $\epsilon =10\mathrm{\%}$, when applied to $\mathrm{\Delta }$ (orange) and $\kappa$ (blue) compared to the ideal designs (dashed). (a) correlation length ${\mathrm{\ell }}_c=\mathrm{\ell }/50$ and (f) ${\mathrm{\ell }}_c=2\mathrm{\ell }$. Other figures give crosstalk $\mu$ versus correlation length with $\epsilon =1\mathrm{\%}$. All couplers have length $\mathrm{\ell }=1$, and, without noise, they are ideal and have no crosstalk. Noise is applied to $\kappa$ only (blue), $\mathrm{\Delta }$ only (orange) and both (yellow) for $N=100$ realisations and averaged. Error bars indicate the standard error of the mean. Crosses on the right are the result as ${\mathrm{\ell }}_c\to \mathrm{\infty }$, calulated using Equation (47)(47) $\begin{array}{rl}& \mu (L_c\to \mathrm{\infty })\\ & =\frac{1}{4\pi {\sigma }^2}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}exp\left(-\frac{{(\delta \kappa )}^2+{(\delta \mathrm{\Delta })}^2}{2{\sigma }^2}\right)\mu (\mathrm{\Delta }+\delta \mathrm{\Delta },\kappa +\delta \kappa )d(\delta \kappa )d(\delta \mathrm{\Delta }),\end{array}$(47) . The noise floor is approximately ${10}^{-8}$. Left column: STA designs using Lewis-Riesenfeld invariants with third-order polynomial for ${\gamma }_{LR}$ and ${\beta }_{LR}$ [21], with different boundary conditions (Equation (24)(24) ${\left.\frac{d{\beta }_{LR}}{dz}\right|}_{z=0}={\left.-\frac{d{\beta }_{LR}}{dz}\right|}_{z=\mathrm{\ell }}=\frac{C}{\mathrm{\ell }},$(24) ): (b) $C=\pi /2$; (c) $C=\pi$; and (d) $C=3\pi /2$. (e) $\mathrm{\Delta }$-robust coupler [23], designed to be robust to variations in $\mathrm{\Delta }$. Right column: conventional couplers: (g) directional coupler; polynomial designs (h) $P_{-2}$; (i) $P_0$; and (j) $P_2$

Figure 9 (b)–(e) and (g)–(f) give the results of our noise analysis for eight different couplers. It shows the crosstalk for noise on $\kappa$ alone (blue), $\mathrm{\Delta }$ alone (orange) and on both $\kappa$ and $\mathrm{\Delta }$ (yellow) by numerically integrating Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) with the modified parameters. For each coupler, we use the same noise profiles, but scaled appropriately to match Equation (46)(46) $\sigma =\frac{\epsilon }{\mathrm{\ell }}{\int }_0^{\mathrm{\ell }}\mathrm{\Gamma }(z)dz.$(46) . This avoids random variations between couplers. We generate $N=100$ realisations of the random process and show the average, as well as the standard error of the mean as the standard deviation of the sample divided by $N$. When ${\mathrm{\ell }}_c\gg \mathrm{\ell }$, the error becomes a Gaussian distributed constant offset over the length of the device, and the effects of the crosstalk can then be computed directly, without averaging, via

(47) $\begin{array}{rl}& \mu (L_c\to \mathrm{\infty })\\ & =\frac{1}{4\pi {\sigma }^2}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}exp\left(-\frac{{(\delta \kappa )}^2+{(\delta \mathrm{\Delta })}^2}{2{\sigma }^2}\right)\mu (\mathrm{\Delta }+\delta \mathrm{\Delta },\kappa +\delta \kappa )d(\delta \kappa )d(\delta \mathrm{\Delta }),\end{array}$(47)

for the case of noise on both $\mathrm{\Delta }$ and $\kappa$. Here $\mu (\mathrm{\Delta },\kappa )$ is the crosstalk from numerically integrating Equation (2)(2) $\frac{d}{dz}|G⟩=-i\left(\begin{array}{cc}\mathrm{\Delta }& \kappa \\ \kappa & -\mathrm{\Delta }\end{array}\right)|G⟩\equiv -iH|G⟩,$(2) . Similarly, a single integral is evaluated for the case of noise on one Hamiltonian element at a time. The results of such calculations are indicated by crosses in Figure 9.

We select a wide range of coupler designs. The STA couplers selected are taken to be of the invariant approach, and are third-order polynomial designs of ${\gamma }_{LR}$ and ${\beta }_{LR}$ [21]. Boundary conditions Equations (22)(22) ${\left.\frac{d{\gamma }_{LR}}{dz}\right|}_{z=0}={\left.\frac{d{\gamma }_{LR}}{dz}\right|}_{z=\mathrm{\ell }}.$(22) –(25(25) ${\beta }_{LR}(0)={\beta }_{LR}(\mathrm{\ell })=-\pi /2.$(25) ) are imposed, with coefficient $C=\pi /2,\pi ,3\pi /2$ (Figure 9(b)–(d)). Recall that this is equivalent to both the counter-diabatic and transitionless tracking methods, hence demonstrating the behaviour of typical STA couplers. The fourth STA coupler is that of Tseng’s, which has been optimized to be insensitive to $\mathrm{\Delta }$ noise [23] (Figure 9(e)) using a perturbation analysis. We compare the performance of these couplers to other designs: a standard directional coupler (Figure 9(g)), and those designed by polynomials $\{P_{-2},P_0,P_2\}$ (Figure 9(h)–(j)). All couplers are designed to produce perfect power transfer at $\mathrm{\ell }=1$ and without noise, but scaling $\mathrm{\Gamma }$ produces the same results at different device lengths. Unsurprisingly, perhaps, the effect of noise on both $\kappa$ and $\mathrm{\Delta }$ is approximately the superposition of the crosstalk with noise on each individually.

Overall, couplers designed by STA methods (left column) tend to have lower crosstalk at longer correlation lengths, while conventional couplers (right column) perform better at shorter correlation lengths. However, their performance at short correlation lengths does not differ markedly from the other designs. Fabrication errors tend to have short correlation lengths relative to the length of typical couplers. Poulton et al. report tens of nanometers [69]. At such scales, there is no significant gain over a simple directional coupler, which can be further optimized for robustness [70]. The primary motivation behind adiabatic couplers is that they are inherently more robust to variations in device parameters, but it appears that this advantage is lost upon the application of STA–the device in Figure 9(b) has poor performance at all correlation lengths with noise applied to both $\kappa$ and $\mathrm{\Delta }$.

6. Discussion and conclusions

The efficient coupling of light between two waveguides is an important function in photonics, and many approaches have been brought to bear on it. While in this Review we have discussed adiabatic couplers, they are related to directional couplers. We have distinguished two broad measures to judge the response. The first of these is the position of the first zero of the crosstalk, indicating the shortest device that can provide perfect coupling. However, such devices would be expected to be sensitive to fabrication errors and other non-idealities. The second measure is the decay rate of the crosstalk envelope–if this envelope decays rapidly with length then the crosstalk might be expected to be robust against perturbations. These two criteria seem to be mutually exclusive and in fact we can distinguish a continuum of couplers for which the directional coupler and high-order polynomial couplers (with $n>0$) are the extreme cases. These polynomial couplers with $n>0$ are characterized by a high degree of smoothness of the parameters at the two edges of the device and rapidly varying parameters in the centre. Such devices have rapidly decreasing envelopes scaling with ${\mathrm{\ell }}^{-2(n+1)}$, but the first zero requires a fairly long device. In contrast, polynomial couplers with negative $n$ have more severe discontinuities at the edges of the device, but vary more slowly in the centre. They have a slowly decreasing crosstalk envelope, but the first zero in the crosstalk appears at short device lengths. The directional coupler is an extreme example of this–its crosstalk envelope is periodic and does not decrease at all, but the first zero occurs as early as $\mathrm{\ell }=\pi /(2\kappa )$. The FAQUAD approach seems to be a compromise between these two. In a recent development, Wu et al. [71] describe a generalization of the concept of FAQUAD, though this intriguing proposal needs more thorough investigation.

The discussion in the previous paragraph does not explicitly include the STA-based design procedures. In our experience, the practical outcomes of these procedures are coupler designs which have zero crosstalk at a particular, designated length. However, such designs do not address the envelope of the crosstalk, nor is the length where the crosstalk vanishes necessarily particularly short. The advantage, then, perhaps lies in other aspects of the design? One of the advantages of STA-based designs is the improvement in bandwidth, as discussed in Section 4.2.4. The sensitivity to fabrication errors is discussed below.

Ideally, one would model the sensitivity of the design to the most likely type of fabrication errors. Such questions have been considered at a general level in previous papers [23,59], but these were carried over from quantum mechanics, where a physical parameter, for example the frequency offset of a laser, is directly related to the elements of the Hamiltonian. In contrast, in photonics this link is more indirect–an error in the width of a waveguide, for example, affects both $\mathrm{\Delta }$ and $\kappa$. It changes $\mathrm{\Delta }$ since it changes the propagation constant of the waveguide mode; it changes $\kappa$ since the evanescent field of that mode changes, and hence the overlap with the waveguide. In addition to this, these errors depend on the particular fabrication method that is used, and it can even depend on the particular piece of equipment. We therefore limited our analysis to investigating the noise model in Section 5 and considered the sensitivity to noise with various correlation lengths. Figure 9 shows the results of a set of such calculations, with a STA-designed coupler on the left and other couplers on the right. Our conclusion is that while neither of the two classes of couplers is obviously superior to the other over the entire range of correlation lengths, which covers four orders of magnitude, some of the devices perform excellently in more limited contexts. In particular, Tseng’s design in Figure 9(d), which is designed to be insensitive to noise in $\mathrm{\Delta }$ [23], outperforms all other devices for $\mathrm{\Delta }$-noise with long correlation lengths. Thus, this design would lead to excellent devices for which this is the predominant fabrication error. However, as discussed above, the relationship between the fabrication errors and device parameters is complicated. The challenge, then, would be to generalize these very powerful methods to other frequency regimes and to correlated errors in $\mathrm{\Delta }$ and $\kappa$.

In conclusion, while the STA-techniques provide a beautiful and consistent mathematical framework, based on the analysis we have carried out, the case for its use in the design of waveguide couplers is not clear. We surmise that this is because, even though the description of such devices can be cast in a Hamiltonian framework, the nature of likely errors is indirect and complicated and may be specific to particular fabrication methods. Considering conventional design methods, the polynomial approach that we earlier introduced [47], augmented here with negative orders, provides a systematic way to balance the need for a rapidly decreasing envelope ($n$ large and positive) with the need for zero crosstalk for a short device, with FAQUAD devices in between as a compromise. These couplers are characterized by having constant $\mathrm{\Gamma }$; generalization to couplers with varying $\mathrm{\Gamma }$ can be considered, but there seems to be no overwhelming advantage of one over the other. In practice, then, the design chosen is likely to depend on the particular device requirements. As a final comment we note that while we only considered two-mode couplers, STA-based techniques have been applied to a large number of photonic devices including multimode devices [39,44]. Although the conclusion reached here may not apply to all of these, we speculate that the crosstalk envelope, that we introduced here as a performance measure, carries over to other devices since they are a consequence of modal interference effects.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

This work was supported by an Australian Research Council Discovery Early Career Researcher Award (DE200101041) (A.Tu.)

Appendix A1.

Equivalence Proofs

A.1. Invariant selection given a transitionless tracking Hamiltonian

In this section, we verify claims in [21], explicitly proving that the proposed invariant Equation (41) satisfies the invariance condition (11) for any system designed with the transitionless tracking method. First, consider the commutator in Equation (11). Since $H_0$ and $I$ are diagonal in the $\left|n(z)\right.⟩$ basis, they commute, i.e. $[I,H]=[I,H_{\mathrm{T}\mathrm{T}}(z)]$. To determine this, we consider the derivative of the resolution of the identity

(A1) $\frac{\mathrm{\partial }}{\mathrm{\partial }z\text{mathbbm}}1=\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left(\sum _n\left|n(z)\right.⟩\left.⟨n(z)\right|\right)=\left(\sum _n\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)\left.⟨n(z)\right|\right)+\left(\sum _n\left|n(z)\right.⟩\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left.⟨n(z)\right|\right)\right)=0.$(A1)

Returning to $[I,H_{\mathrm{T}\mathrm{T}}(z)]$ we find

$i[I,H_{\mathrm{T}\mathrm{T}}(z)]=$

$=\left(\sum _n\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)\left.⟨n(z)\right|\right)\left(\sum _n\left|n(z)\right.⟩{\lambda }_n\left.⟨n(z)\right|\right)-\left(\sum _n\left|n(z)\right.⟩{\lambda }_n\left.⟨n(z)\right|\right)\left(\sum _n\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)\left.⟨n(z)\right|\right)$

(A2) $=\left(\sum _n\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)\left.⟨n(z)\right|\right)\left(\sum _n\left|n(z)\right.⟩{\lambda }_n\left.⟨n(z)\right|\right)+\left(\sum _n\left|n(z)\right.⟩{\lambda }_n\left.⟨n(z)\right|\right)\left(\sum _n\left|n(z)\right.⟩\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left.⟨n(z)\right|\right)\right),$(A2)

Where we used a previous result to change from the derivative of $\left|n(z)\right.⟩$ to $\left.⟨n(z)\right|$. Simplifying gives

(A3) $[I,H_{\mathrm{T}\mathrm{T}}(z)]=-i\left(\sum _n{\lambda }_n\left(\left|n(z)\right.⟩\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left.⟨n(z)\right|\right)+\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)\left.⟨n(z)\right|\right)\right).$(A3)

The first term of Equation (11)(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11) is given by

(A4) $\frac{\mathrm{\partial }}{\mathrm{\partial }z}I=\left(\sum _n{\lambda }_n\left(\left|n(z)\right.⟩\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left.⟨n(z)\right|\right)+\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n(z)\right.⟩\right)\left.⟨n(z)\right|\right)\right).$(A4)

Hence, the invariance condition Equation (11)(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11) is satisfied for any choice of coupler designed by the transitionless tracking method.

Appendix A2

A.2. Invariant selection given a unitary transformed transitionless tracking Hamiltonian

In this section, we follow a similar line of reasoning to Appendix A1, showing explicitly that Equation (43) is the photonics equivalent of the invariant proposed by [21]. We consider a general unitary transformation $V$, since this property turns out to be quite general. We first consider the derivative of the identity

(A5) $\frac{\mathrm{\partial }}{\mathrm{\partial }z}\text{mathbbm}1=\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left(V^{†}(z)V(z)\right)=\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}^{†}(z)\right)V(z)+V^{†}(z)\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}V(z)\right)=0.$(A5)

We now drop the explicit $z$ dependence for clarity, recalling that only ${\lambda }_n$ is independent of $z$. Consider the commutator, which will have three terms once Equation (37) is used

(A6) $[I_V,H_V]=[I_V,V{H}_0{V}^{†}]+[I_V,V{H}_{TT}{V}^{†}]-i\left[I_V,V\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}^{†}\right)\right].$(A6)

The first of these is zero, since they both share the $V\left|n(z)\right.⟩$ basis. Computing the second term with the definition of $H_{TT}$ from Equation (37) gives

(A7) $\begin{array}{rl}& -i[I_V,V{H}_{TT}{V}^{†}]=\\ & =\left(\sum _{n}V\left|n\right.⟩{\lambda }_n⟨n|V^{†}\right)\left(\sum _{n}V\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n\right.⟩\right)⟨n|V^{†}\right)-\left(\sum _{n}V\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n\right.⟩\right)⟨n|V^{†}\right)\\ & \left(\sum _{n}V\left|n\right.⟩{\lambda }_n⟨n|V^{†}\right)\\ & =-\left(\sum _{n}V\left|n\right.⟩{\lambda }_n⟨n|V^{†}\right)\left(\sum _{n}V\left|n\right.⟩\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}⟨n|\right)V^{†}\right)-\left(\sum _{n}V\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n\right.⟩\right)⟨n|V^{†}\right)\\ & \left(\sum _{n}V\left|n\right.⟩{\lambda }_n⟨n|V^{†}\right)\\ & =-\sum _n\left(V\left|n\right.⟩{\lambda }_n\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}⟨n|\right)V^{†}\right)-\sum _n\left(V\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n\right.⟩\right){\lambda }_n⟨n|V^{†}\right),\end{array}$(A7)

where we have again used the derivative of the identity (A1) to change the derivative from $\left|n(z)\right.⟩$ to $\left.⟨n(z)\right|$. Similarly, for the third commutator term

(A8) $\begin{array}{rl}& \left[I_V,V\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}^{†}\right)\right]=\left(\sum _{n}V\left|n\right.⟩{\lambda }_n⟨n|V^{†}\right)\left(V\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}^{†}\right)\right)\\ & -\left(V\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}^{†}\right)\right)\left(\sum _{n}V\left|n\right.⟩{\lambda }_n⟨n|V^{†}\right)\\ & =\left(\sum _{n}V\left|n\right.⟩{\lambda }_n⟨n|\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}^{†}\right)\right)-\left(\sum _{n}V\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}^{†}\right)V\left|n\right.⟩{\lambda }_n⟨n|V^{†}\right)\\ & =\left(\sum _{n}V\left|n\right.⟩{\lambda }_n⟨n|\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}^{†}\right)\right)+\left(\sum _n\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}V\right)\left|n\right.⟩{\lambda }_n⟨n|V^{†}\right),\end{array}$(A8)

where we have used the result (A5). The remaining term in Equation (11)(11) $\frac{dI}{dz}\equiv \frac{\mathrm{\partial }I}{\mathrm{\partial }z}+\frac{1}{i}[I(z),H(z)]=0,$(11) is given by

(A9) $\frac{\mathrm{\partial }{I}_V}{\mathrm{\partial }z}=\sum _n{\lambda }_n\left(\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}V\right)\left|n\right.⟩⟨n|V^{†}+V\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}\left|n\right.⟩\right)⟨n|V^{†}+V\left|n\right.⟩\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}⟨n|\right)V^{†}+V\left|n\right.⟩⟨n|\left(\frac{\mathrm{\partial }}{\mathrm{\partial }z}{V}^{†}\right)\right).$(A9)

We observe that all terms in Equation (11) cancel, proving that Equation (43) is an invariant for a system under the unitary transformed dynamics described in Section 3.4.