Solution of Cross-Kerr Interaction Combined with Parametric Amplification

We present a full operator approach to treatment of the cross-Kerr interaction combined with parametric amplification. It is shown that this problem can be exactly integrated using the method of higher-order operators. While the initial basis is infinite-dimensional, an orthogonal transformation can reduce the problem exactly into a six-dimensional basis which can be integrated conveniently.


S2 Silence of Multiplicative Operators
The first issue to notice in treatment of the multiplicative noise is the dependence of the multiplying operators to the noise terms. These operators also are determined from lower order Langevin equations in which similar noise terms are fed in. Iteratively going back to the lowest order determines that these multiplying operators appear as an infinite series such aŝ a(t) =â 0 (t) + a(t) + a 1 (t)â in (t) + a 2 (t)â 2 in (t) + . . . , whereâ 0 (t) = exp [(i[M] − [Γ])t]â 0 (0) is the decaying operator term of the homogeneous solution to the system of Langevin equations, which decays to zero, thus taking no part in the steady state solution, and is excluded from contributing to the noise spectral density. Furthermore, a(t) = a(t)1 is the noiseless or silent part of the operatorâ(t) driven by the external classical pump field, which can in principle be determined from solving the Langevin equations with the stochastic terms dropped, while only keeping the drive terms as input. The above term ultimately gets multiplied to another Gaussian white noise term, such asâ in (t) again. Such a multiplicative noise term asâ(t)â in (t) will have an expansion given bŷ and so on. It is not difficult to see that as long as a lower-order Gaussian noise term is present, the higher-order terms will have negligible contribution to non-zero absolute (and not detuned) frequencies in the ultimate noise spectral density. This is discussed in §S5.
In order to establish this, we may define the second-order noise corresponding to the squared processĉ(t) = 1 2â 2 (t), which clearly has a decay rate of 2κ. The corresponding stochastic noise process iŝ The stochastic processĉ in (t) is no longer Gaussian white althoughâ in (t) is a Gaussian white stochastic process by assumption with a symmetrized auto-correlation â † (t)â(τ) S = 1 2 δ (t − τ). The symmetrized autocorrelation of this higher-order stochastic process in light of the Isserlis-Wick theorem [S1] is thus given by

S1/S4
The corresponding spectral density of this noise process, being its Fourier transform, simply causes a Dirac delta at zero frequency [S1]. Similarly, all higher-power noise processes will have no contribution to the non-zero frequency of the noise spectral density. As a result, the multiplicative noise (S2) can be effectively truncated aŝ without causing any error in the non-zero frequencies of the resulting noise spectral density. A more general treatment of the second-order noise processes with Gaussian resonances is discussed elsewhere [S1].

S3 Noise Spectral Density
Following the general approach to construction of the scattering matrix based on the input-output formalism [S1-S3], one may easily show that where S(w) is the 6 × 6 scattering matrix given by Expansion of the output operator array giveŝ whereÂ in,l (w) are multiplicative noise terms such as a l (w)â in,l (w), whereâ in,l stand for white Guassian White stochastic processesâ in ,b in , and their conjugatesb † in ,b † in , and also a l (w) are the corresponding Fourier-transformed silent multiplicative terms. This can be correspondingly shown to lead to the noise spectral densities with the understanding that the terms corresponding to conjugate noise operators are grouped together under the absolute value.
Here, S A l A l (w) are the symmetrized noise spectral densities of the Gaussian White processesâ in,l . The spectral densities of these processes are typically constants asn l + 1 2 withn l being thermal occupation number of bosons. For an optical bosonic bath, one may conveniently setn l = 0, while for phononsn l can be estimated from Bose-Einstein distribution [S4]. Furthermore, the symbol * represents convolution in the frequency domain.
Neither of these is the directly measurable spectrum, but it is rather the noise spectral density of ladder operatorb for photons, which can be measured. These necessitates a way to recover the information through what is calculable by the method of higher-order operators. The symmetrized noise spectral density ofd = 1 2b 2 is by definition given in terms of the Fourier transform of the corresponding symmetrized auto-correlation function, which is where the last expression is found by application of the Isserlis-Wick theorem and b † (τ)b † (τ) S = 0. By noting the definition of Fourier and inverse Fourier transforms, we get where 1 2 is substrated and added to account for the half a quanta of white noise which is lost in the symmetrization, and if not removed will cause appearance of a non-physical Dirac delta under the square root. While S DD (w) is found from simple scattering matrix calculations, all it takes now to find the measurable quantity S BB (w) is to take an inverse Fourier transform, followed by a square root and another Fourier transform. Similarly, one we may now defineS BB (w) = 1 2 [S BB (w) + S BB (−w)] as the symmetrized noise spectrum. The equation (S11) is the main key to recover the expected results from the higher-order operator algebra.

S4 Time-evolution of Operators
It is easy to obtain the explicit solution to the truncated system of Langevin equations (5) where exp(·) represents the matrix exponentiation, and

S5 Non-negative Integer Powers of Noise
It is straightforward to see that any term involving a non-negative integer power of a noise such asα j j in (t) has identically zero contribution to the measured noise spectral density. In order to show this, let us the noise assume the normal autocorrelation In the limit of ζ → +∞ this will settle back to the expected Dirac delta's function δ (t − τ). The autocorrelation of the measurable optical field is connected to the operatorâ(t), which by means of the Isserlis-Wick theorem becomes The corresponding noise spectral density in frequency domain, where w is the absolute optical frequency (and not the detuning referenced to a certain non-zero resonance frequency), is given by In the limit of ζ → +∞ with j > 1 the above expression is identically zero, and hence meeting the claim. It is equally straightbackward to show that for any white Gaussian noise such asâ in satisfying â † in (t)â in (τ) S = 1 2 δ (t − τ), the higher-power noise processesα j in (t) = κ 1− j 2â j in (t) contribute only to the zero frequency of the noise spectral density. To show this, we assume which again in the limit of ζ → +∞ reproduces the Dirac's delta δ (t − τ). Then the Isserlis-Wick theorem for such a Gaussian noise process could be exactly used to write Taking the Fourier transform from both sides gives the resulting noise spectral density which in the limit of κ = ζ → +∞ yields an upper bound to a constant number of quanta √ j/2 j , being less than 1 2 for j > 1. This maximum bound to the background number of added noise quanta due to higher-power noise rapidly decays to zero with increasing j.
These limits are physically meaningful as long as cavity linewidth is much larger than the pump laser linewidth, which is quite accurately met in practice. So, when no squeezing is taking place and cavity resonances exhibit noise spectra corresponding to a much larger number of quanta than 1 2 , it should be safe to ignore the effect of square noise terms and higher powers. Ultimately, a numerical integration carried out on a nonlinear differential equation with exaggerated noise input amplitude could very visibly distinguish the zero contribution of the higher-power noise terms, quite expectedly, confirming the general above conclusions.

S6 Multiplicative Noise
It is possible to make simple estimates for b(w) to be used instead ofb in (30), where a convolution such as (S9) would have been needed instead. This can be done by setting up the Langevin equations for operators {b,b † }, which are coupled. Then these have to be Fourier-transformed and diagonalized to findb(w) explicitly, and the expectation of this expression b(w) = b (w) can be now used as [S5] S DD (w) = 1 γ 2 |b(w) * [S 11 (w) + S 13 (w)] + b * (w) * [S 12 (w) + S 13 (w)]| 2 S YY (w).
This process is a lot more complicated and typically can be simplified by direct utilization of (30). Nevertheless, the corresponding Langevin equations are Replacingn withn, taking the expectation and some simplification gives [W] = 1 2 Here, the factor 1 2 is included to take care of normalization with respect to 2Ω rather than Ω. Functions b(w) and b * (w) from solution of (S21) can be plugged in the convolutions of (S19). Note the cancellation of γ as it explicitly appears in (S19) and (S21).