Implementation of bright six-partite entanglement by coupled intracavity sum frequency generation

Quantum entanglement is the foundation of quantum physics and the key resource of quantum information processing [1]. Over the course of the past few decades, continuous-variable (CV) entanglement, especially multipartite CV entanglement, has drawn much attention in its potential applications in quantum communication and computation, which include quantum-controlled dense coding [2, 3], quantum telecloning [4, 5] and quantum teleportation network [6]. Multipartite CV entanglement fields with different frequencies are more useful in connecting different physical systems at the nodes of quantum networks in contrast to those with the same frequency generated by beam splitters [7, 8]. Many papers have reported the schemes which can generate two-color or multi-color CV entangled beams by an optical parametric down-conversion process [9–13], frequency up-conversion process [14–17] and cascaded or concurrent nonlinear processes [18–29].

In addition, the nonlinear optical coupler consists of two parallel nonlinear optical waveguides, which coupled by evanescent overlaps of the guided modes can generate two-color or multi-color CV entanglement. This type of coupling was theoretically and experimentally investigated decades ago [30, 31]. The coupled degenerate parametric down-conversion can generate bipartite [32–34] or two-color quadripartite CV entanglement [35] in traveling-wave configuration or in an optical cavity. Moreover, the coupled non-degenerate parametric down-conversion in an optical cavity can produce bright three-color six-partite entanglement [36]. Olsen et al studied the Kerr nonlinear coupler that produces bipartite CV entangled beams and can be used to demonstrate the Einstein–Podolsky–Rosen (EPR) paradox [37]. Coupled up-conversion processes, unlike coupled parametric down-conversion processes, can also generate bright multipartite CV entanglement. Since second harmonic generation (SHG) and sum frequency generation (SFG) have no oscillation threshold, the generated entangled fields become macroscopically occupied as soon as the low-frequency (LF) pumping fields are turned on. Bache et al theoretically discussed the coupled SHGs and analyzed the intensity correlations as well as the phase-dependent correlation features [38]. Bright bipartite [39] and two-color quadripartite [40] entanglement can be produced by coupled SHGs inside an optical cavity. The related system is all-integrated, which makes it a compact robust source of CV entanglement.

In this present letter, we proposed a scheme for the generation of bright three-color six-partite CV entanglement from the coupled SFG processes in an optical cavity in which periodically poled LiTaO₃ waveguides serve as the gain medium. The periods can be designed flexibly according to the practical application. The quantum correlations among these six output beams and the steady-state regions are analyzed according to the inseparability criterion for multipartite CV entanglement as proposed by van Loock and Furusawa [41].

The proposed scheme, which consists of two parallel adjacent nonlinear waveguides and an optical cavity, is shown in figure 1. The two waveguides are fabricated in a periodically poled quadratic nonlinear optical crystal. Efficient SFG processes take place in each waveguide. These two waveguides are numbered A and B, respectively. M1 and M2 are the coupler mirrors. The linear coupling between these waveguides is achieved by the evanescent wave overlaps of the intracavity guide modes. LF pump fields $a_{1}^{{\rm in}},b_{1}^{{\rm in}} ~$ (with the frequency of ${{\omega}_{1}}$) and $a_{2}^{{\rm in}},~b_{2}^{{\rm in}}$ (with the frequency of ${{\omega}_{2}}$) are incident upon the nonlinear waveguides A and B through the coupler mirror M1. ${{a}_{0}}, ~ {{b}_{0}}$ (with the frequency of ${{\omega}_{0}}$) are sum frequency (SF) fields in waveguides A and B, respectively. The energy conservation condition is ${{\omega}_{1}}+{{\omega}_{2}}={{\omega}_{0}}$. Here, we assume that the two SFG processes are equal and all the modes inside the waveguides are perfectly phase matched using the technique of quasi-phase-matching (QPM). Because the periods can be flexibly designed according to QPM, almost any desired wavelength of the entangled beams can be obtained by modification of the phase-matching conditions and the LF pump wavelengths. For example, we can choose two LF pump field wavelengths as 1560 and 780 nm, respectively. The SF field is centered at 520 nm. 1560 nm is suitable to be transmitted through optical fibers over long distances and 780 nm is one of the absorption lines of the rubidium atom, which is an ideal candidate for the storage in quantum repeater of quantum information.

Zoom In Zoom Out Reset image size

The effective Hamiltonian for this system can be written as

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\widehat{H}}_{{\rm tol}}}={{\widehat{H}}_{{\rm sys}}}+{{\widehat{H}}_{{\rm couple}}}+{{\widehat{H}}_{{\rm bath}}}.\nonumber \end{align} \tag{ 1 }$$

The interaction Hamiltonian for the coupled system is

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \begin{array}{@{}cccccccccccccccccccc@{}} & {{\widehat{H}}_{{\rm sys}}}={\rm i}\hbar \kappa \left(\widehat{a}_{1}^{+}\widehat{a}_{2}^{+}{{\widehat{a}}_{0}}-{{\widehat{a}}_{1}}{{\widehat{a}}_{2}}\widehat{a}_{0}^{+}+\widehat{b}_{1}^{+}\widehat{b}_{2}^{+}{{\widehat{b}}_{0}}-{{\widehat{b}}_{1}}{{\widehat{b}}_{2}}\widehat{b}_{0}^{+} \right)-\hbar [{{\delta}_{1}}(\widehat{a}_{1}^{+}{{\widehat{a}}_{1}} \nonumber \\ & +\,\widehat{b}_{1}^{+}{{\widehat{b}}_{1}})+{{\delta}_{2}}\left(\widehat{a}_{2}^{+}{{\widehat{a}}_{2}}+\widehat{b}_{2}^{+}{{\widehat{b}}_{2}} \right)+{{\delta}_{0}}\left(\widehat{a}_{0}^{+}{{\widehat{a}}_{0}}+\widehat{b}_{0}^{+}{{\widehat{b}}_{0}} \right)] \end{array}\nonumber \end{align} \tag{ 2 }$$

where, the first part of the Hamiltonian stands for the SFG processes, and the second part is the free Hamiltonian of the cavity modes for the LP pump and SF fields. $\kappa$ is the dimensionless nonlinear coupling coefficient of our system, which is proportional to the nonlinear susceptibility and the structure parameters of the periodically poled nonlinear optical waveguides, and taken to be real without loss of generality [9]. ${{\widehat{a}}_{n}}$, and ${{\widehat{b}}_{n}}~(n=0,1,2)$ are bosonic annihilation operators for the mode at frequency ${{\omega}_{n}}~\left(n=0,1,2 \right)$ in waveguides A and B, respectively. ${{\delta}_{n}}={{\omega}_{n}}-\omega _{n}^{{\rm cav}}~(n=0,1,2)$ are the cavity detunings from their respective resonances.

The coupling Hamiltonian is described by

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\widehat{H}}_{{\rm couple}}}=&\,\hbar {{J}_{1}}\left(\widehat{a}_{1}^{+}{{\widehat{b}}_{1}}+\widehat{b}_{1}^{+}{{\widehat{a}}_{1}} \right)+\hbar {{J}_{2}}\left(\widehat{a}_{2}^{+}{{\widehat{b}}_{2}}+\widehat{b}_{2}^{+}{{\widehat{a}}_{2}} \right)\nonumber \\ &+\hbar {{J}_{0}}\left(\widehat{a}_{0}^{+}{{\widehat{b}}_{0}}+\widehat{b}_{0}^{+}{{\widehat{a}}_{0}} \right).\nonumber \end{align} \tag{ 3 }$$

In the above, ${{J}_{n}}~\left(n=0,1,2 \right)$ are the linear coupling parameters at three different frequencies of two nonlinear waveguides. As [38] described, the lower frequency (${{\omega}_{1}}$ and ${{\omega}_{2}}$) coupling parameters ${{J}_{1}}$ and ${{J}_{2}}$ are larger than the higher frequency (${{\omega}_{0}}$) coupling parameter ${{J}_{0}}$.

Finally, the cavity damping can be written as

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\widehat{H}}_{{\rm bath}}}=\hbar \underset{n=0}{\overset{2}{\mathop \sum}}\,[\widehat{\Gamma}_{a}^{n}\widehat{a}_{n}^{+}+\widehat{\Gamma}_{b}^{n}\widehat{b}_{n}^{+}]+{\rm H}.{\rm c}.,\nonumber \end{align} \tag{ 4 }$$

where $\widehat{\Gamma}_{a}^{n}$ and $\widehat{\Gamma}_{b}^{n}~(n=0,1,2)$ are the bath operators in two waveguides for three different frequencies, respectively.

The quantum Langevin equations of motion for the six intracavity fields according to the input–output formalism developed by Collett and Gardiner [42] can be obtained as

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \begin{array}{@{}cccccccccccccccccccc@{}} & \tau \frac{{\rm d}{{\widehat{a}}_{1}}}{{\rm d}t}=(-{{\gamma}_{1}}+{\rm i}{{\Delta}_{1}}){{\widehat{a}}_{1}}+\kappa \widehat{a}_{2}^{+}{{\widehat{a}}_{0}}-{\rm i}{{J}_{1}}{{\widehat{b}}_{1}}+\sqrt{2{{\gamma}_{1}}}\widehat{a}_{1}^{{\rm in}}, \nonumber \\ & \tau \frac{{\rm d}{{\widehat{a}}_{2}}}{{\rm d}t}=(-{{\gamma}_{2}}+{\rm i}{{\Delta}_{2}}){{\widehat{a}}_{2}}+\kappa \widehat{a}_{1}^{+}{{\widehat{a}}_{0}}-{\rm i}{{J}_{2}}{{\widehat{b}}_{2}}+\sqrt{2{{\gamma}_{2}}}\widehat{a}_{2}^{{\rm in}}, \nonumber \\ & \tau \frac{{\rm d}{{\widehat{a}}_{0}}}{{\rm d}t}=(-{{\gamma}_{0}}+{\rm i}{{\Delta}_{0}}){{\widehat{a}}_{0}}-\kappa {{\widehat{a}}_{1}}{{\widehat{a}}_{2}}-{\rm i}{{J}_{0}}{{\widehat{b}}_{0}}+\sqrt{2{{\gamma}_{0}}}\widehat{a}_{0}^{{\rm in}}, \nonumber \\ & \tau \frac{{\rm d}{{\widehat{b}}_{1}}}{{\rm d}t}=(-{{\gamma}_{1}}+{\rm i}{{\Delta}_{1}}){{\widehat{b}}_{1}}+\kappa \widehat{b}_{2}^{+}{{\widehat{b}}_{0}}-{\rm i}{{J}_{1}}{{\widehat{a}}_{1}}+\sqrt{2{{\gamma}_{1}}}\widehat{b}_{1}^{{\rm in}}, \nonumber \\ & \tau \frac{{\rm d}{{\widehat{b}}_{2}}}{{\rm d}t}=\left(-{{\gamma}_{2}}+{\rm i}{{\Delta}_{2}} \right){{\widehat{b}}_{2}}+\kappa \widehat{b}_{1}^{+}{{\widehat{b}}_{0}}-{\rm i}{{J}_{2}}{{\widehat{a}}_{2}}+\sqrt{2{{\gamma}_{2}}}\widehat{b}_{2}^{{\rm in}}, \nonumber \\ & \tau \frac{{\rm d}{{\widehat{b}}_{0}}}{{\rm d}t}=\left(-{{\gamma}_{0}}+{\rm i}{{\Delta}_{0}} \right){{\widehat{b}}_{0}}-\kappa {{\widehat{b}}_{1}}{{\widehat{b}}_{2}}-{\rm i}{{J}_{0}}{{\widehat{a}}_{0}}+\sqrt{2{{\gamma}_{0}}}\widehat{b}_{0}^{{\rm in}} \end{array}\nonumber \end{align} \tag{ 5 }$$

where $\tau$ is the cavity round-tip time and assumed to be the same for all the cavity modes. $\widehat{a}_{n}^{{\rm in}}~(n=0,1,2)$ and $\widehat{b}_{n}^{{\rm in}}$ are the input field operators of the cavity of two waveguides. The dimensionless variables ${{\Delta}_{n}}=~{{\delta}_{n}}\tau ~(n=0,1,2)$, represent the cavity detunings of the SF and two LF pump fields, respectively. ${{\gamma}_{n}}(n=0,1,2)$ are the loss coefficients of the cavity modes, which are related to the amplitude reflection coefficients and the amplitude transmission coefficients.

By setting the left side of equation (5) to be zero and replacing all the operators with their mean values, we can get a set of equations for the classical steady-state values, i.e.

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \begin{array}{@{}cccccccccccccccccccc@{}} & 0=\left(-{{\gamma}_{1}}+{\rm i}{{\Delta}_{1}} \right){{\alpha}_{1}}+\kappa \alpha _{2}^{{\rm *}}{{\alpha}_{0}}-{\rm i}{{J}_{1}}{{\beta}_{1}}+\sqrt{2{{\gamma}_{1}}}\alpha _{1}^{{\rm in}}, \nonumber \\ & 0=\left(-{{\gamma}_{2}}+{\rm i}{{\Delta}_{2}} \right){{\alpha}_{2}}+\kappa \alpha _{1}^{{\rm *}}{{\alpha}_{0}}-{\rm i}{{J}_{2}}{{\beta}_{2}}+\sqrt{2{{\gamma}_{2}}}\alpha _{2}^{{\rm in}}, \nonumber \\ & 0=\left(-{{\gamma}_{0}}+{\rm i}{{\Delta}_{0}} \right){{\alpha}_{0}}-\kappa {{\alpha}_{1}}{{\alpha}_{2}}-{\rm i}{{J}_{0}}{{\beta}_{0}}+\sqrt{2{{\gamma}_{0}}}\alpha _{0}^{{\rm in}}, \nonumber \\ & 0=\left(-{{\gamma}_{1}}+{\rm i}{{\Delta}_{1}} \right){{\beta}_{1}}+\kappa \beta _{2}^{{\rm *}}{{\beta}_{0}}-{\rm i}{{J}_{1}}{{\alpha}_{1}}+\sqrt{2{{\gamma}_{1}}}\beta _{1}^{{\rm in}}, \nonumber \\ & 0=\left(-{{\gamma}_{2}}+{\rm i}{{\Delta}_{2}} \right){{\beta}_{2}}+\kappa \beta _{1}^{{\rm *}}{{\beta}_{0}}-{\rm i}{{J}_{2}}{{\alpha}_{2}}+\sqrt{2{{\gamma}_{2}}}\beta _{2}^{{\rm in}}, \nonumber \\ & 0=\left(-{{\gamma}_{0}}+{\rm i}{{\Delta}_{0}} \right){{\beta}_{0}}-\kappa {{\beta}_{1}}{{\beta}_{2}}-{\rm i}{{J}_{0}}{{\alpha}_{0}}+\sqrt{2{{\gamma}_{0}}}\beta _{0}^{{\rm in}}, \end{array}\nonumber \end{align} \tag{ 6 }$$

where ${{\alpha}_{n}}(0,1,2)$ and ${{\beta}_{n}}$ are the steady-state amplitudes of the intracavity modes ${{\widehat{a}}_{n}}$ and ${{\widehat{b}}_{n}}$, respectively. For simplicity, we will make some assumptions: two waveguides are equal and the initial phases of two coherent LF pump fields are to be zero. Therefore, we set $\sqrt{2{{\gamma}_{1}}}\alpha _{1}^{{\rm in}}=\sqrt{2{{\gamma}_{1}}}\beta _{1}^{{\rm in}}={{\varepsilon}_{1}}$ and $\sqrt{2{{\gamma}_{2}}}\alpha _{2}^{{\rm in}}=\sqrt{2{{\gamma}_{2}}}\beta _{2}^{{\rm in}}={{\varepsilon}_{2}}$ and no input of the SF fields, i.e. $\sqrt{2{{\gamma}_{0}}}\alpha _{0}^{{\rm in}}=\sqrt{2{{\gamma}_{0}}}\alpha \beta _{0}^{{\rm in}}=0$. Also, the steady-state solutions in two waveguides are equal, i.e. ${{\alpha}_{1}}={{\beta}_{1}}={{A}_{1}}{{{\rm e}}^{{\rm i}{{\phi}_{1}}}}$, ${{\alpha}_{2}}={{\beta}_{2}}={{A}_{2}}{{{\rm e}}^{{\rm i}{{\phi}_{2}}}}$ and ${{\alpha}_{0}}={{\beta}_{0}}={{A}_{0}}{{{\rm e}}^{{\rm i}{{\phi}_{0}}}}$.

In the following, we will discuss this system in the case of detuning the cavity by an appropriate amount from three frequencies so as to simplify the theoretical analysis and actually improve some quantum correlations, i.e. $~{{\Delta}_{1}}={{J}_{1}}$, ${{\Delta}_{2}}={{J}_{2}}$ and ${{\Delta}_{0}}={{J}_{0}}$. Then, we will solve the steady-state equations in two cases, the first symmetric case, where the frequencies of two LF pump fields are close to each other, i.e. ${{\omega}_{1}}$ is close to ${{\omega}_{2}}$, and we can assume the system parameters of these modes in two waveguides as ${{\gamma}_{1}}={{\gamma}_{2}}=\gamma$, ${{\varepsilon}_{1}}={{\varepsilon}_{1}}=\varepsilon$, $\alpha ={{\alpha}_{1}}={{\alpha}_{2}}={{\beta}_{1}}={{\beta}_{2}}=A{{{\rm e}}^{{\rm i}\phi}}$, ${{\Delta}_{1}}={{\Delta}_{2}}=\Delta$, ${{J}_{1}}={{J}_{2}}=J$ and $\Delta =J$. Therefore, under this assumption, we can get the steady-state solutions for the symmetry case as

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \begin{array}{@{}cccccccccccccccccccc@{}} & \alpha =\frac{\varepsilon}{\gamma -\kappa {{\alpha}_{0}}} \nonumber \\ & {{\alpha}_{0}}=\frac{1}{6}\left(\frac{4\gamma}{\kappa}+\frac{{{2}^{\frac{2}{3}}}A}{{{\gamma}_{0}}{{\kappa}^{2}}}+\frac{{{16}^{1/3}}{{\gamma}^{2}}{{\gamma}_{0}}}{A} \right) \nonumber \\ & {{\phi}_{0}}-{{\phi}_{1}}-{{\phi}_{2}}=\pi , ~ ~ {{\phi}_{1}}={{\phi}_{2}}=0 \end{array}\nonumber \end{align} \tag{ 7 }$$

with $A=\sqrt[3]{-2{{\gamma}^{3}}\gamma _{0}^{3}{{\kappa}^{3}}-27{{\varepsilon}^{2}}{{\kappa}^{5}}\gamma _{0}^{2}+\sqrt{27{{\varepsilon}^{2}}{{\kappa}^{8}}\gamma _{3}^{4}(4{{\gamma}^{3}}{{\gamma}_{0}}+27{{\varepsilon}^{2}}{{\kappa}^{2}})}}$.

The other asymmetric case, i.e. ${{\gamma}_{1}}\ne {{\gamma}_{2}}$, ${{\varepsilon}_{1}}\ne {{\varepsilon}_{2}}$, ${{\alpha}_{1}}\ne {{\alpha}_{2}}$, ${{\beta}_{1}}\ne {{\beta}_{2}}$, ${{\Delta}_{1}}\ne {{\Delta}_{2}}$ and ${{J}_{1}}\ne {{J}_{2}}$. Therefore, the steady-state solutions for the asymmetry case can be written as

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \begin{array}{@{}cccccccccccccccccccc@{}} & {{A}_{0}}=\frac{\kappa {{A}_{1}}{{A}_{2}}}{{{\gamma}_{0}}} \nonumber \\ & {{A}_{2}}=\frac{{{\varepsilon}_{2}}}{{{\gamma}_{2}}+{{\kappa}^{2}}A_{1}^{2}/{{\gamma}_{0}}} \nonumber \\ & ~ {{\phi}_{0}}-{{\phi}_{1}}-{{\phi}_{2}}=\pi , ~ {{\phi}_{1}}={{\phi}_{2}}=0 \end{array}\nonumber \end{align} \tag{ 8 }$$

with ${{A}_{1}}$ being the solution of the following high-order equation:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \left({{\varepsilon}_{1}}-{{\gamma}_{1}}{{A}_{1}} \right){{({{\gamma}_{2}}{{\gamma}_{0}}+{{\kappa}^{2}}A_{1}^{2})}^{2}}-{{\gamma}_{0}}{{\kappa}^{2}}\varepsilon _{2}^{2}{{A}_{1}}=0.\nonumber \end{align} \tag{ 9 }$$

Because equation (9) is a five-order equation, we cannot get its analytical solutions of steady-state equations in the asymmetric case.

In order to analyze the spectral correlations of the output modes, we may linearize the system variables into their classical steady-state solutions and small quantum fluctuations, i.e. ${{\widehat{a}}_{n}}={{\alpha}_{n}}+\delta {{\widehat{a}}_{n}}$ and ${{\widehat{b}}_{n}}={{\beta}_{n}}+\delta {{\widehat{b}}_{n}} ~ (n=0,1,2)$. Here, we define $\delta {{\widehat{X}}_{an}}=\delta {{\widehat{a}}_{n}}+\delta \widehat{a}_{n}^{+}$ and $\delta {{\widehat{X}}_{bn}}=\delta {{\widehat{b}}_{n}}+\delta \widehat{b}_{n}^{+}~(n=0,1,2)$, as the fluctuations of the amplitude quadrature and $\delta {{\widehat{Y}}_{an}}=-{\rm i}(\delta {{\widehat{a}}_{n}}-\delta \widehat{a}_{n}^{+})$ and $\delta {{\widehat{Y}}_{bn}}=-{\rm i}(\delta {{\widehat{b}}_{n}}-\delta \widehat{b}_{n}^{+})$ stand for the fluctuations of the phase quadrature of each mode. Substituting these expressions into equation (5) and neglecting second- and higher-order terms of the fluctuations, we can obtain the linearized equation for the fluctuations as

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \tau \frac{{\rm d}D}{{\rm d}t}=MD+B{{D}^{{\rm in}}}\nonumber \end{align} \tag{ 10 }$$

with

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle D={{[\delta {{\widehat{X}}_{a1}},\delta {{\widehat{Y}}_{a1}},\delta {{\widehat{X}}_{a2}},\delta {{\widehat{Y}}_{a2}},\delta {{\widehat{X}}_{a0}},\delta {{\widehat{Y}}_{a0}},\delta {{\widehat{X}}_{b1}},\delta {{\widehat{Y}}_{b1}},\delta {{\widehat{X}}_{b2}},\delta {{\widehat{Y}}_{b2}},\delta {{\widehat{X}}_{b0}},\delta {{\widehat{Y}}_{b0}}]}^{T}},\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{D}^{{\rm in}}}={{[\delta \widehat{X}_{a1}^{{\rm in}},\delta \widehat{Y}_{a1}^{{\rm in}},\delta \widehat{X}_{a2}^{{\rm in}},\delta \widehat{Y}_{a2}^{{\rm in}},\delta \widehat{X}_{a0}^{{\rm in}},\delta \widehat{Y}_{a0}^{{\rm in}},\delta \widehat{X}_{b1}^{{\rm in}},\delta \widehat{Y}_{b1}^{{\rm in}},\delta \widehat{X}_{b2}^{{\rm in}},\delta \widehat{Y}_{b2}^{{\rm in}},\delta \widehat{X}_{b0}^{{\rm in}},\delta \widehat{Y}_{b0}^{{\rm in}}]}^{T}},\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle M=\left(\begin{array}{@{}cccccccccccccccccccc@{}} {{M}_{1}} & {{M}_{2}} \nonumber \\ {{M}_{2}} & {{M}_{1}} \end{array}\!\! \right),\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{M}_{1}}=\left(\begin{array}{@{}cccccccccccccccccccc@{}} -{{\gamma}_{1}} & -{{\Delta}_{1}} & {{m}_{0}} & {{n}_{0}} & {{m}_{2}} & {{n}_{2}} \nonumber \\ {{\Delta}_{1}} & -{{\gamma}_{1}} & {{n}_{0}} & -{{m}_{0}} & -{{n}_{2}} & {{m}_{2}} \nonumber \\ {{m}_{0}} & {{n}_{0}} & -{{\gamma}_{2}} & -{{\Delta}_{2}} & {{m}_{1}} & {{n}_{1}} \nonumber \\ {{n}_{0}} & -{{m}_{0}} & {{\Delta}_{2}} & -{{\gamma}_{2}} & -{{n}_{1}} & {{m}_{1}} \nonumber \\ -{{m}_{2}} & {{n}_{2}} & -{{m}_{1}} & {{n}_{1}} & -{{\gamma}_{0}} & -{{\Delta}_{0}} \nonumber \\ -{{n}_{2}} & -{{m}_{2}} & -{{n}_{1}} & -{{m}_{1}} & {{\Delta}_{0}} & -{{\gamma}_{0}} \end{array}\!\! \right),\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{M}_{2}}=\left(\begin{array}{@{}cccccccccccccccccccc@{}} 0 & {{J}_{1}} & 0 & 0 & 0 & 0 \nonumber \\ -{{J}_{1}} & 0 & 0 & 0 & 0 & 0 \nonumber \\ 0 & 0 & 0 & {{J}_{2}} & 0 & 0 \nonumber \\ 0 & 0 & -{{J}_{2}} & 0 & 0 & 0 \nonumber \\ 0 & 0 & 0 & 0 & 0 & {{J}_{0}} \nonumber \\ 0 & 0 & 0 & 0 & -{{J}_{0}} & 0 \end{array} \right),\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle B=\left(\begin{array}{@{}cccccccccccccccccccc@{}} {{B}_{1}} & 0 \nonumber \\ 0 & {{B}_{1}} \end{array} \right),\quad {{B}_{1}}=\left(\begin{array}{@{}cccccccccccccccccccc@{}} \sqrt{2{{\gamma}_{1}}} & 0 & 0 & 0 & 0 & 0 \nonumber \\ 0 & \sqrt{2{{\gamma}_{1}}} & 0 & 0 & 0 & 0 \nonumber \\ 0 & 0 & \sqrt{2{{\gamma}_{2}}} & 0 & 0 & 0 \nonumber \\ 0 & 0 & 0 & \sqrt{2{{\gamma}_{2}}} & 0 & 0 \nonumber \\ 0 & 0 & 0 & 0 & \sqrt{2{{\gamma}_{0}}} & 0 \nonumber \\ 0 & 0 & 0 & 0 & 0 & \sqrt{2{{\gamma}_{0}}} \end{array} \right)\nonumber \end{align}$$

where ${{m}_{0}}=\kappa {{A}_{0}}{\rm cos}{{\phi}_{0}}$, ${{n}_{0}}=\kappa {{A}_{0}}{\rm sin}{{\phi}_{0}}$, ${{m}_{1}}=\kappa {{A}_{1}}{\rm cos}{{\phi}_{1}}$, ${{n}_{1}}=\kappa {{A}_{1}}{\rm sin}{{\phi}_{1}}$, ${{m}_{2}}=\kappa {{A}_{2}}{\rm cos}{{\phi}_{2}}$ and ${{n}_{2}}=\kappa {{A}_{2}}{\rm sin}{{\phi}_{2}}$.

As long as none of the eigenvalues of matrix M develop a positive real part, the fluctuations will not tend to increase exponentially, and the linearized fluctuation analysis is valid. In this present letter, we only study the region where the linearization analysis is valid. Because the analytical eigenvalues of matrix M cannot be obtained, we give the numerical results of the real parts of the eigenvalues. For the symmetry case, the real parts of the eigenvalues of matrix M are schematically shown in figure 2 versus different system parameters $\gamma$, ${{\gamma}_{0}}$, $\varepsilon$, $\kappa$ and $J$. In figures 2(a) and (b), we plot the real parts of the eigenvalues versus the loss coefficients $\gamma$ and ${{\gamma}_{0}}$ with the system parameters chosen as $\varepsilon =400$, $\kappa =0.01$, $J=10\gamma$ and ${{J}_{0}}=\gamma /4$. Figure 2(c) depicts the dimensionless nonlinear coupling coefficients $\kappa$ with $\gamma =2$ and ${{\gamma}_{0}}=10$. Figures 2(d) and (e) depict the real parts of the eigenvalues versus the symmetric pump power ε and the coupling parameter $J$ of two LF fields, respectively.

Zoom In Zoom Out Reset image size

For the asymmetry case, one can numerically obtain the real parts of the eigenvalues of drift matrix M versus different system parameters ${{\gamma}_{1}}$, ${{\gamma}_{2}}$, ${{\gamma}_{0}}$, ${{\varepsilon}_{1}}$, ${{\varepsilon}_{2}}$, $\kappa$, ${{J}_{1}}$ and ${{J}_{2}}$ as shown in figure 3. The parameters used in figures 3(a)–(c) are ${{\varepsilon}_{1}}=860$, ${{\varepsilon}_{2}}=1320$, $\kappa =0.01$, ${{J}_{1}}=10{{\gamma}_{1}}$, ${{J}_{2}}=4{{\gamma}_{1}}$ and ${{J}_{0}}={{\gamma}_{1}}/4$. In figure 3(d), the parameters are chosen as ${{\gamma}_{1}}=2$, ${{\gamma}_{2}}=14$ and ${{\gamma}_{0}}=2$, while the other parameters remain unchanged. In figures 3(e)–(h) $\kappa =0.01$ and the rest of the parameters are the same as in figure 3(d).

Zoom In Zoom Out Reset image size

From figures 2 and 3, one can see that the values of the solid lines, which represent the real parts of the eigenvalues of drift matrix M are all negative in a wide range of all given system parameters, which indicate clearly that the system is stable and the linearization analysis is valid in this region. In this present letter, we only study the region where the linearization analysis is valid.

By using the Fourier transformation and the boundary condition ${{D}^{{\rm out}}}=\sqrt{2{{\gamma}_{n}}}D-{{D}^{{\rm in}}}~\left(n=0,1,2 \right)$ on the coupler mirror [42], the output fluctuations in the frequency space can be obtained as

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{D}^{{\rm out}}}=\left[B{{\left({\rm i}\omega \tau I-M \right)}^{-1}}B-I \right]{{D}^{{\rm in}}}\nonumber \end{align} \tag{ 11 }$$

where $I$ is the unit matrix.

In the following, we will investigate the entanglement characteristics of six output beams of this system by employing the sufficient inseparability criterion for multipartite CV entanglement, which is given by van Loock and Furusawa [41]. In order to sufficiently verify that the six beams of this system are fully inseparable and possess the property of genuine six-partite entanglement, the following five inequalities should be violated at the same time:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \begin{array}{@{}cccccccccccccccccccc@{}} & {{S}_{a1,a2}}={{\delta}^{2}}\left({{\widehat{X}}_{a1}}+{{\widehat{X}}_{a2}} \right)+{{\delta}^{2}}\left({{\widehat{Y}}_{a1}}-{{\widehat{Y}}_{a2}}+{{g}_{3}}{{\widehat{Y}}_{a0}}+{{g}_{4}}{{\widehat{Y}}_{b1}}+{{g}_{5}}{{\widehat{Y}}_{b2}}+{{g}_{6}}{{\widehat{Y}}_{b0}} \right)\geqslant 4 \nonumber \\ & {{S}_{a1,a0}}={{\delta}^{2}}\left({{\widehat{X}}_{a1}}-{{\widehat{X}}_{a0}} \right)+{{\delta}^{2}}\left({{\widehat{Y}}_{a1}}+{{g}_{2}}{{\widehat{Y}}_{a2}}+{{\widehat{Y}}_{a0}}+{{g}_{4}}{{\widehat{Y}}_{b1}}+{{g}_{5}}{{\widehat{Y}}_{b2}}+{{g}_{6}}{{\widehat{Y}}_{b0}} \right)\geqslant 4 \nonumber \\ & {{S}_{b1,b2}}={{\delta}^{2}}\left({{\widehat{X}}_{b1}}+{{\widehat{X}}_{b2}} \right)+{{\delta}^{2}}\left({{g}_{1}}{{\widehat{Y}}_{a1}}+{{g}_{2}}{{\widehat{Y}}_{a2}}+{{g}_{3}}{{\widehat{Y}}_{a0}}+{{\widehat{Y}}_{b1}}-{{\widehat{Y}}_{b2}}+{{g}_{6}}{{\widehat{Y}}_{b0}} \right)\geqslant 4 \nonumber \\ & {{S}_{b1,b0}}={{\delta}^{2}}\left({{\widehat{X}}_{b1}}-{{\widehat{X}}_{b0}} \right)+{{\delta}^{2}}\left({{g}_{1}}{{\widehat{Y}}_{a1}}+{{g}_{2}}{{\widehat{Y}}_{a2}}+{{g}_{3}}{{\widehat{Y}}_{a0}}+{{\widehat{Y}}_{b1}}+{{g}_{5}}{{\widehat{Y}}_{b2}}+{{\widehat{Y}}_{b0}} \right)\geqslant 4 \nonumber \\ & {{S}_{a1,b1}}={{\delta}^{2}}\left({{\widehat{X}}_{a1}}-{{\widehat{X}}_{b1}} \right)+{{\delta}^{2}}\left({{\widehat{Y}}_{a1}}+{{g}_{2}}{{\widehat{Y}}_{a2}}+{{g}_{3}}{{\widehat{Y}}_{a0}}+{{\widehat{Y}}_{b1}}+{{g}_{5}}{{\widehat{Y}}_{b2}}+{{g}_{6}}{{\widehat{Y}}_{b0}} \right)\geqslant 4. \end{array}\nonumber \end{align} \tag{ 12 }$$

In order to study the correlation between two SF fields of two waveguides, we also take the next equation into account:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{S}_{a0,b0}}={{\delta}^{2}}\left({{\widehat{X}}_{a0}}-{{\widehat{X}}_{b0}} \right)+{{\delta}^{2}}\left({{g}_{1}}{{\widehat{Y}}_{a1}}+{{g}_{2}}{{\widehat{Y}}_{a2}}+{{\widehat{Y}}_{a0}}+{{g}_{4}}{{\widehat{Y}}_{b1}}+{{g}_{5}}{{\widehat{Y}}_{b2}}+{{\widehat{Y}}_{b0}} \right)\geqslant 4\nonumber \end{align}$$

where ${{g}_{s}}~(s=1,2,3,4,5,6)$ are scaling factors, which can be freely adjusted to minimize the correlation spectra functions.

We will first numerically discuss the quantum correlation spectra in the symmetry case versus different system parameters. Due to the symmetry of two waveguides, we have ${{S}_{a1,a2}}={{S}_{b1,b2}}$ and ${{S}_{a1,a0}}={{S}_{b1,b0}}$. In figure 4, we plot the quantum correlation spectra versus the analysis frequency $\Omega =\omega \tau /{{\gamma}_{1}}$ with $\varepsilon =400$, $\gamma =2$, ${{\gamma}_{0}}=6$, $J=10\gamma$, ${{J}_{0}}=\gamma /4$ and $\kappa =0.01$, respectively. As shown in this figure, these five inequalities are violated, so we can conclude that genuine six-partite entanglement can be achieved. In addition, the correlations between two LF pump fields are better than those of any other fields, and the best six-partite entanglement is obtained at about $\Omega =0.62$.

Zoom In Zoom Out Reset image size

The dependence of quantum correlation spectra on LF pump power $\varepsilon$ is shown in figure 5, where $\gamma =2$, ${{\gamma}_{0}}=6$, $J=10\gamma$, ${{J}_{0}}=\gamma /4$, $\kappa =0.01$ and $\Omega =0.65$, respectively. In this figure, one can see that the values of these five curves are all less than four when $\varepsilon$ is greater than 170, which demonstrates that these six modes are entangled. For the given set of parameters, the best six-partite entanglement is achieved around $\varepsilon =405$. Also, the correlations between the LF pump fields is much better than those of other fields, and the correlations between the generated SF fields are first increased and then decreased as the LF pump fields increase. This characteristic is the same as the uncoupled SFG.

Zoom In Zoom Out Reset image size

Figure 6 depicts the quantum correlation spectra versus the linear coupling parameter of two LF pump fields $J$ with $\varepsilon =400$, $\gamma =2$, ${{\gamma}_{0}}=6$, $J=10\gamma$, ${{J}_{0}}=\gamma /4$, $\kappa =0.01$ and $\Omega =0.65$, respectively. As shown in this figure, when the value of $J$ is very small, i.e. the linear interaction between the two waveguides is small, there is only quantum correlation among the modes of each waveguide, but no correlation between the modes from different waveguides and these six modes are not entangled with each other. With increasing the value of $J$ and when $J>0.4$, the values of all the correlation curves are less than four. Continuously increasing the value of $J$, one can find that almost all the values of the correlation spectra do not change, which indicates that these six modes are entangled with each other as long as $J>5$.

Zoom In Zoom Out Reset image size

Next, we will discuss the asymmetry case, i.e. ${{\gamma}_{1}}\ne {{\gamma}_{2}}$, ${{\varepsilon}_{1}}\ne {{\varepsilon}_{2}}$, ${{\alpha}_{1}}\ne {{\alpha}_{2}}$, ${{\beta}_{1}}\ne {{\beta}_{2}}$, ${{\Delta}_{1}}\ne {{\Delta}_{2}}$ and ${{J}_{1}}\ne {{J}_{2}}$. We plot the quantum correlation spectra versus the analysis frequency $\Omega =\omega \tau /{{\gamma}_{1}}$ in figure 7 with ${{\varepsilon}_{1}}=860$, ${{\varepsilon}_{2}}=1320$, ${{\gamma}_{1}}=1$, ${{\gamma}_{2}}=14$, ${{\gamma}_{0}}=8$, ${{J}_{1}}={{\Delta}_{1}}=10{{\gamma}_{1}}$, ${{J}_{2}}={{\Delta}_{2}}=4{{\gamma}_{1}}$, ${{J}_{0}}={{\Delta}_{0}}={{\gamma}_{1}}/4$ and $~\kappa =0.01$, respectively. It is obvious that the values of these six quantum correlation spectra are all less than four, and violate six inequalities in equation (12). Therefore, these six output fields are genuine six-partite entanglement in the asymmetrical case. The best correlation can be obtained at $\Omega =0$.

Zoom In Zoom Out Reset image size

Figure 8 shows the result of the quantum correlation spectra versus the ratio of two LF pump power $\sigma ={{\varepsilon}_{2}}/{{\varepsilon}_{1}}$ with ${{\varepsilon}_{1}}=860$, $\Omega =0$ and the other parameters remain unchanged as in figure 7. As shown in this figure, when $\sigma <0.07$ there is no entanglement between any two fields. This might be because the power of pump field ${{\varepsilon}_{2}}$ is much lower than that of ${{\varepsilon}_{1}}$, and the quantum property of pump field ${{\varepsilon}_{1}}$ is not present in such a case. By increasing the power of pump field ${{\varepsilon}_{2}}$, the nonlinear conversion efficiency is increased, while pump field ${{\varepsilon}_{1}}$ becomes weaker, and two LF pump fields begin to entangle with each other. With further increasing the power, the quantum properties of these six fields are present and bright three-color six-partite CV entanglement appears when $\sigma >1.08$. Moreover, the correlation between LF pump fields is always better than that of pump ${{\varepsilon}_{1}}$ and the SF field. The correlation of the fields in the same waveguide is better than that of fields in different waveguides. Compared with other cases, we achieve the largest degree of entanglement when $\sigma =2.93$.

Zoom In Zoom Out Reset image size

Figure 9 depicts the relationship between the quantum correlation spectra with the linear coupling parameter ${{J}_{1}}$. We set ${{\varepsilon}_{1}}=860$ and ${{\varepsilon}_{2}}=2.93{{\varepsilon}_{1}}$ and other parameters remain unchanged. From this figure, we can see that the result is similar to figure 6; i.e. when the value of ${{J}_{1}}$ is very small, there is no entanglement among the output fields of different waveguides, and it only exists when the quantum correlation is among the fields of the same waveguide. When ${{J}_{1}}$ is increased to ${{J}_{1}}=0.63$, the values of all the correlation curves are less than four and, six-partite entanglement appears. And as long as ${{J}_{1}}>2$ almost all the values of the correlation spectra do not change.

Zoom In Zoom Out Reset image size

In conclusion, we proposed a theoretical scheme for producing bright three-color six-partite CV entanglement based on coupled SFG processes. LF pump fields exit the optical cavity from one side, while spatially separated and SF beams exit the optical cavity from the other side. We theoretically analyzed the quantum correlation spectra of six output fields in symmetric and asymmetric cases and found that the entanglement characteristics are influenced by different system parameters. This integrated scheme may have applications in the field of quantum cryptography and quantum communication networks.

This research was supported by the NNSFC (grant nos. 11504175, 11021403), the Special Funds of NNSFC (grant no. 51245010), the NSF of Jiangsu Province (grant no. BK2012463), the open research fund program of Nanjing University of Information Science and Technology (grant no. 16KF071) and the College Students Practice Innovation Training Program of Nuist (grant nos. 201710300261 and 1171130001001).