format: text

text:
# Disagreement between correlations of quantum mechanics and stochastic electrodynamics in the damped parametric oscillator

## I Introduction

Local hidden variable theories have been extensively compared to quantum mechanics over the last seventy or so years . Most comparisons between the two have investigated whether or not quantum mechanics is equivalent to a local hidden variable theory. Much evidence indicates that it is not. Many results in quantum mechanics have been found that are incompatible with all local hidden variable theories . Most of these results have involved idealized undamped systems. However, all experimental systems encounter damping. Thus, it is interesting (and more realistic) to compare quantum mechanics and local hidden variable theories in damped systems . This paper compares quantum mechanics and one local hidden variable theory (SED) in such a system.

One of the earliest works comparing local hidden variable theories to quantum mechanics was Bell’s theorem . It demonstrates that quantum mechanics is incompatible with all local hidden variable theories at a statistical level. It does so by deriving an upper bound on a function of two particle correlations for all local hidden variable theories, which quantum mechanics exceeds. Extensions of it have been formulated for large angular momentum and particle number systems . These extensions demonstrate nonclassical behavior in a regime usually regarded as being purely classical. Greenberger, Horne and Zeilinger (GHZ)  have also extended Bell’s work, differentiating quantum mechanics from all local hidden variable theories for single, as opposed to ensemble, measurements. The three particle GHZ theorem has an “all or nothing” quality and distinguishes between local hidden variable theories and quantum mechanics in a single experimental run, once three basic correlations are established.

Comparisons between quantum mechanics and local hidden variable theories have also been made using continuous variables (which are discretized in formulating the comparison), such as quadrature phase amplitudes  and there is currently much interest in this area. For quadrature phase amplitude measurements, these comparisons can have detector efficiencies of in excess of 99% . They also tend to relate more strongly to Einstein, Podolsky and Rosen’s original EPR paradox  than earlier discrete variable ones. Indeed, the EPR paradox has been experimentally demonstrated using quadrature phase amplitudes . Additionally, quantum teleportation has been achieved using quadrature phase amplitudes  further demonstrating the utility of continuous variables.

One commonly used local hidden variable theory is stochastic electrodynamics (SED) . Some authors have proposed it as an alternative to quantum mechanics . Furthermore, a semiclassical approach equivalent to it  is also commonly used in parametric oscillator calculations . SED consists of adding Gaussian white noise to classical electrodynamics. It is equivalent to truncating third order derivative terms in the quantum mechanical Moyal equation, a commonly used approximation . Such terms are often negligible and thus SED reproduces many results of quantum mechanics . However, it cannot violate Bell inequalities for quadrature phase amplitude measurements, and is thus distinct from quantum mechanics . Various authors have explicitly shown differences between SED and quantum mechanics . In particular, it has been shown that the two theories predict different transient third order correlations for the undamped nondegenerate parametric oscillator . It has also been shown that they predict different macroscopic quadrature phase amplitude correlations in the damped nondegenerate parametric oscillator in the steady state .

In general, differences between quantum mechanics and local hidden variable theories are reduced or eliminated by damping . Furthermore, damping is a significant element of many realistic systems. It is thus important to consider its effects on differences between quantum mechanics and local hidden variable theories such as SED. However, all but a few of the comparisons between the quantum mechanics and local hidden variable theories referenced above have involved undamped systems. They are thus idealized in this respect. In contrast, damping is included in the calculatons in this paper. It is included to consider a theoretical model which is as realistic as possible and also to determine the sensitivity of differences between quantum mechanics and SED to its presence.

This paper extends a previous comparison between quantum mechanics and SED in the nondegenerate parametric oscillator . In particular, it contrasts both intracavity and external moments of the two theories’ in the same system with damping included. Expressions from both theories are compared for the intracavity moment $`\mathrm{\Delta }X_1(\tau )\mathrm{\Delta }(X_2(\tau )\mathrm{\Delta }X_3(\tau )`$, where $`\mathrm{\Delta }X_i(\tau )=X_i(\tau )X_i(\tau )`$, for $`i=1,2,3`$, $`X_i(\tau )`$ is quadrature phase amplitude, the subscripts represent different radiation modes and $`\tau `$ is a scaled time variable. A comparison is also made for an analogous external moment. Both analytic iterative and numerical techniques are used to calculate moments. The results produced by these techniques show that the intracavity and external moments differ greatly between the two theories. In particular, the analytic method shows that the moments of quantum mechanics are cubic in the system’s nonlinear coupling constant to leading order whilst those of SED are linear. The two theories are compared over a range of nonlinear coupling constant, damping and average initial pump photon number values. The results of these comparisons show a number of qualitative trends. Most importantly, quantum mechanics and SED differ in the situations considered with the largest particle number and damping to nonlinear coupling ratios, although the differences are reduced in relative size.

Stochastic techniques are used to obtain results both for quantum mechanics and SED. The positive-P coherent state representation  is used to calculate quantum mechanical predictions. It is particularly well suited to the calculation of quantum dynamics in damped quantum optical systems when nonclassical behavior is present. It is able to handle arbitrarily large photon numbers. It converges quickly (in the sense of sampling error) when systems’ dimensionless nonlinearities are relatively small, as is the case with nonlinear optical experiments. By contrast, the method used for SED calculations corresponds to commonly used approaches in quantum optics, where the field is treated as a semiclassical object surrounded by (classical) vacuum fluctuations. Both methods are used to generate analytic predictions and are also numerically simulated.

## II Quantum mechanics

This paper considers an idealized nondegenerate parametric oscillator, resonant at three frequencies $`\omega _1\mathrm{and}\omega _2`$ (signal and idler frequencies) and $`\omega _3=\omega _1+\omega _2`$ (pump frequency). It contains a nonlinear medium that couples the modes and converts higher energy pump photons into lower energy signal and idler ones. The system’s interaction Hamiltonian, including linear losses, is given by

$$\widehat{H}=i\mathrm{}G(\widehat{a}_1^{}\widehat{a}_2^{}\widehat{a}_3\widehat{a}_1\widehat{a}_2\widehat{a}_3^{})+\underset{i=1}{\overset{3}{}}\widehat{\mathrm{\Gamma }}_i\widehat{a}_i^{}+\widehat{\mathrm{\Gamma }}_i^{}\widehat{a}_i,$$
(1)
where $`\widehat{a}_i^{}\mathrm{and}\widehat{a}_i`$ are creation and annihilation operators for oscillator modes, $`\widehat{\mathrm{\Gamma }}_i^{}\mathrm{and}\widehat{\mathrm{\Gamma }}_i`$ are environment mode operators and G is a nonlinear interaction strength constant. Initially, the system has a coherent state in the pump mode, and vacuum states in the signal and idler modes.

A number of quasiprobability representations exist to describe quantum states, the most famous being the Glauber-Sudarshan representation. It is produced by decomposing quantum density operators using a diagonal coherent state basis. Thus,

$$\widehat{\rho }=𝑑\alpha ^2P(\alpha ,\alpha ^{})|\alpha \alpha |,$$
(2)
where $`\widehat{\rho }`$ is a density operator and $`P(\alpha ,\alpha ^{})`$ is the Glauber-Sudarshan representation. The Glauber-Sudarshan representation can be negative and is hence not a strict probability density function. A more recent representation is the positive-P representation  which is an actual probability density function over an off-diagonal coherent state basis. It further differs from the Glauber-Sudarshan representation by using a phase space of doubled dimension. The positive-P variables $`\{\alpha _i,\alpha _i^+\}`$, where $`i`$ is a positive integer, are analogous to complex field amplitudes, with $`\alpha _i`$ and $`\alpha _i^+`$ describing a particular radiation mode. However, $`\{\alpha _i\}`$ and $`\{\alpha _i^+\}`$ are independent and hence $`\alpha _i(\alpha _i^+)^{}`$, though their averages are complex conjugate and thus $`\alpha _i=\alpha _i^+^{}`$. Variable averages are equal to normally ordered quantum averages once the substitutions $`\alpha _i\widehat{a}_i`$ and $`\alpha _i^+\widehat{a}_i^{}`$ are made. For example, $`\alpha _1\alpha _1^+=\widehat{a}_1^{}\widehat{a}_1_\rho `$, where $`\widehat{O}_\rho `$ denotes $`Tr(\widehat{\rho }\widehat{O})`$, as usual in quantum mechanics.

Stochastic equations of motion for positive-P variables for the damped nondegenerate parametric oscillator are, in terms of $`\tau `$ (time scaled by $`\mathrm{\Gamma }`$, a typical damping constant with units of inverse time),

$`{\displaystyle \frac{d\alpha _1}{d\tau }}`$ $`=`$ $`\gamma _1\alpha _1+g\alpha _2^+\alpha _3+\sqrt{g\alpha _3}\xi _1`$ (3)
$`{\displaystyle \frac{d\alpha _1^+}{d\tau }}`$ $`=`$ $`\gamma _1\alpha _1^++g\alpha _2\alpha _3^++\sqrt{g\alpha _3^+}\xi _1^+`$ (4)
$`{\displaystyle \frac{d\alpha _2}{d\tau }}`$ $`=`$ $`\gamma _2\alpha _2+g\alpha _1^+\alpha _3+\sqrt{g\alpha _3}\xi _2`$ (5)
$`{\displaystyle \frac{d\alpha _2^+}{d\tau }}`$ $`=`$ $`\gamma _2\alpha _2^++g\alpha _1\alpha _3^++\sqrt{g\alpha _3^+}\xi _2^+`$ (6)
$`{\displaystyle \frac{d\alpha _3}{d\tau }}`$ $`=`$ $`\gamma _3\alpha _3g\alpha _1\alpha _2`$ (7)
$`{\displaystyle \frac{d\alpha _3^+}{d\tau }}`$ $`=`$ $`\gamma _3\alpha _3^+g\alpha _1^+\alpha _2^+.`$ (8)
Here $`\xi _1,\xi _2,\xi _1^+\mathrm{and}\xi _2^+`$ are complex Gaussian white noises with the following correlations

$`\xi _i(\tau _1)\xi _j(\tau _2)`$ $`=`$ $`\delta _{3i,j}\delta (\tau _1\tau _2)`$ (9)
$`\xi _i^+(\tau _1)\xi _j^+(\tau _2)`$ $`=`$ $`\delta _{3i,j}\delta (\tau _1\tau _2)`$ (10)
$`\xi _i^+(\tau _1)\xi _j(\tau _2)`$ $`=`$ $`0,`$ (11)
where $`i,j=1,2`$. In Eq. (3), $`\gamma _i=\mathrm{\Gamma }_i/\mathrm{\Gamma }`$, where $`\mathrm{\Gamma }_i`$ is a damping constant for mode $`i`$ with units of inverse time, $`g=G/\mathrm{\Gamma }`$ and $`\tau =\mathrm{\Gamma }t`$. It is assumed that G, $`\mathrm{\Gamma }_i`$ and $`\mathrm{\Gamma }`$ are real. Initial conditions are $`\alpha _1(0)=0,\alpha _2(0)=0\mathrm{and}\alpha _3(0)=ϵ`$. It is noted that Eq. (3) is only valid when boundary terms in phase space can be neglected. These are asymptotically small in the limit of short times or large damping ratios .

Eq. 3 is solved using an analytic iterative method. This method treats damping terms exactly, and noise and nonlinear terms iteratively. It involves, firstly, rewriting the equations forming Eq. (3) as $`\dot{\alpha }_i=\gamma _i\alpha _i+f_i(\{\alpha _j,\alpha _j^+\},\tau )`$ or $`\dot{\alpha }_i^+=\gamma _i\alpha _i^++f_i^+(\{\alpha _j,\alpha _j^+\},\tau )`$, where $`i,j=1,2,3`$. Successively higher order approximations for $`\{\alpha _i(\tau ),\alpha _i^+(\tau )\}`$ are then found using increasingly better approximations for $`f_i`$ and $`f_i^+`$. Thus, $`(m+1)^{th}`$ order terms are given by

$`\alpha _i^{(m+1)}(\tau )`$ $`=`$ $`\alpha _i^{(0)}(\tau )+{\displaystyle _{\tau _1=0}^{\tau _1=\tau }}𝑑\tau _1\mathrm{exp}[\gamma _i(\tau _1\tau )]f_i(\{\alpha _j^{(m)},\alpha _j^{+(m)}\},\tau _1)`$ (12)
$`\alpha _i^{+(m+1)}(\tau )`$ $`=`$ $`\alpha _i^{+(0)}(\tau )+{\displaystyle _{\tau _1=0}^{\tau _1=\tau }}𝑑\tau _1\mathrm{exp}[\gamma _i(\tau _1\tau )]f_i^+(\{\alpha _j^{(m)},\alpha _j^{+(m)}\},\tau _1),`$ (13)
where $`\alpha _k^{(0)}(\tau )=\alpha _k(\tau =0)\mathrm{exp}(\gamma _k\tau )`$ and $`\alpha _k^{(0)+}(\tau )=\alpha _k^{(0)}(\tau )^{}`$ where $`k=1,2,3`$. For example,

$$\alpha _1^{(m+1)}(\tau )=\alpha _1^{(0)}(\tau )+_{\tau _1=0}^{\tau _1=\tau }𝑑\tau _1\mathrm{exp}[\gamma _1(\tau _1\tau )]\left(g\alpha _2^{+(m)}(\tau _1)\alpha _3^{(m)}(\tau _1)+\sqrt{g\alpha _3^{(m)}(\tau _1)}\xi _1(\tau _1)\right)$$
(14)
and first order approximations are

$`\alpha _i^{(1)}(\tau )`$ $`=`$ $`{\displaystyle _{\tau _i=0}^{\tau _i=\tau }}𝑑\tau _i\mathrm{exp}[\gamma _i(\tau _i\tau )]\sqrt{gϵ}\mathrm{exp}({\displaystyle \frac{\gamma _i\tau _i}{2}})\xi _i(\tau _1)`$ (15)
$`\alpha _3^{(1)}(\tau )`$ $`=`$ $`ϵ\mathrm{exp}(\gamma _3\tau )`$ (16)
$`\alpha _i^{+(1)}(\tau )`$ $`=`$ $`{\displaystyle _{\tau _i=0}^{\tau _i=\tau }}𝑑\tau _i\mathrm{exp}[\gamma _i(\tau _i\tau )]\sqrt{gϵ^{}}\mathrm{exp}({\displaystyle \frac{\gamma _i\tau _i}{2}})\xi _i^+(\tau _1)`$ (17)
$`\alpha _3^{+(1)}(\tau )`$ $`=`$ $`ϵ^{}\mathrm{exp}(\gamma _3\tau ),`$ (18)
where $`i=1,2`$.

## III Stochastic diagrams

The iterative method of the previous section can be used, in conjunction with stochastic diagrams, to readily produce analytic approximations for the intracavity moments of quantum mechanics considered in this paper. Stochastic diagrams  are schematic representations of the combinatoric parts of an iterative process. They clearly lay out all terms produced by different orders of iteration. Fundamental stochastic diagrams appear as one of three classes. Those associated with initial conditions appear as straight lines, those with noise terms as straight lines with a cross at their end and those with nonlinear terms as straight lines containing a fork, as shown in Fig. 1(a)-(c). Higher order iterative terms are represented by stochastic diagrams using either combinations of the three basic classes. For example, one of the iterative terms in $`\alpha _1^{(2)}(\tau )`$ is

$$_{\tau _1=0}^{\tau _1=\tau }𝑑\tau _1\mathrm{exp}[\gamma _1(\tau _1\tau )]g\alpha _3^{(0)}(\tau _1)_{\tau _2=0}^{\tau _2=\tau _1}𝑑\tau _2\mathrm{exp}[\gamma _2(\tau _2\tau _1)]\sqrt{gϵ^{}}\mathrm{exp}(\frac{\gamma _2\tau _2}{2})\xi _2^+(\tau _2).$$

It combines all three basic classes and is represented by the stochastic diagram in Fig. 1(d). All iterative terms can be represented by stochastic diagrams.

Stochastic diagrams can also be used to determine the orders of iterative terms. In particular, they can be used to determine the orders of such terms in the system’s nonlinear coupling constant g. This paper focuses on the order of terms in this constant. For quantum mechanics, initial value iterative terms are $`O(g^0)`$, noise iterative terms $`O(g^{1/2})`$ and nonlinear iterative terms $`O(g)`$. Hence, lines in stochastic diagrams count as order zero, crosses as order 1/2 and vertices as order 1. A term’s order is simply found by considering its stochastic diagram and adding a half to its order for every cross and one for every vertex. For example, the term represented in Fig. 1(d) has one vertex and one cross and thus is $`O(g^{\frac{3}{2}})`$. A notation that denotes the order in g of a term by a superscript \[n\] is used in this section.

Stochastic diagrams are now used to determine the intracavity moments of quantum mechanics considered in this paper. Consider all eight moments of the form $`\mathrm{\Delta }𝒜_1(\tau )\mathrm{\Delta }𝒜_2(\tau )\mathrm{\Delta }𝒜_3(\tau )`$, where $`\mathrm{\Delta }𝒜_\mathrm{i}(\tau )=𝒜_i(\tau )𝒜_i(\tau )`$ and $`𝒜_i(\tau )`$ is either $`\widehat{a}_i`$ or $`\widehat{a}_i^{}`$. These are equal to the positive-P variable moments which replace $`\widehat{a}_i`$ and $`\widehat{a}_i^{}`$ by $`\alpha _i(\tau )`$ and $`\alpha _i^+(\tau )`$ respectively. Now, consider the equations that constitute Eq. (3). Their forms do not change when they are expressed in terms of $`\alpha _i(\tau ),\alpha _i^+(\tau ),\alpha _3(\tau )`$ and $`\alpha _3^+(\tau )`$, where $`i=1,2`$. From this, it follows that $`𝒜_i(\tau )=𝒜_i(\tau )`$ and hence $`\alpha _i(\tau )=\alpha _i^+(\tau )=0`$, where again $`i=1,2`$. Thus, $`\mathrm{\Delta }𝒜_1(\tau )\mathrm{\Delta }𝒜_2(\tau )\mathrm{\Delta }𝒜_3(\tau )`$ can be simplified to $`𝒜_1(\tau )𝒜_2(\tau )\mathrm{\Delta }𝒜_3(\tau )`$.

An approximate expression for $`𝒜_1(\tau )𝒜_2(\tau )\mathrm{\Delta }𝒜_3(\tau )`$ is now obtained using the iterative method in Section II (and stochastic diagrams). This method can be used to produce power series expressions in $`g`$ for the positive-P variables. These expressions can then be used to generate power series expressions in $`g`$ for the moments of the form $`𝒜_1(\tau )𝒜_2(\tau )\mathrm{\Delta }𝒜_3(\tau )`$. As $`g1`$ in realistic systems, these power series expressions can be approximated by their lowest order nonzero terms.

Figs 2 (a) and (b) show the stochastic diagrams required to determine the moments of the form $`𝒜_1(\tau )𝒜_2(\tau )\mathrm{\Delta }𝒜_3(\tau )`$. Naively, it might be thought that the lowest order nonzero terms from $`𝒜_1(\tau ),𝒜_2(\tau )\mathrm{and}\mathrm{\Delta }𝒜_3(\tau )`$ simply need to be multiplied together and the average of the subsequent product determined to calculate the lowest order nonzero term in $`𝒜_1(\tau )𝒜_2(\tau )\mathrm{\Delta }𝒜_3(\tau )`$. This is not always true. Sometimes, $`𝒜_1(\tau ),𝒜_2(\tau )\mathrm{and}\mathrm{\Delta }𝒜_3(\tau )`$ are not necessarily zero and yet $`𝒜_1(\tau )𝒜_2(\tau )\mathrm{\Delta }𝒜_3(\tau )`$ is zero. For example, the lowest order nonzero terms for the positive-P variables in $`\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3^+(\tau )`$ are

$$\alpha _i^{[1/2]}(\tau )=_{\tau _i=0}^{\tau _i=\tau }𝑑\tau _i\mathrm{exp}[\gamma _i(\tau _i\tau )]\sqrt{gϵ}\mathrm{exp}(\frac{\gamma _i\tau _i}{2})\xi _i(\tau _i),$$
(19)
where $`i=1,2`$ and

$`\mathrm{\Delta }\alpha _3^{+[2]}(\tau )`$ $`=`$ $`{\displaystyle _{\tau _3=0}^{\tau _3=\tau }}𝑑\tau _3\mathrm{exp}[\gamma _3(\tau _3\tau )]g{\displaystyle _{\tau _4=0}^{\tau _4=\tau _3}}𝑑\tau _4\mathrm{exp}[\gamma _1(\tau _4\tau _3)]\sqrt{gϵ^{}}\mathrm{exp}({\displaystyle \frac{\gamma _1\tau _4}{2}})\xi _1^+(\tau _4)`$ (23)
$`{\displaystyle _{\tau _5=0}^{\tau _5=\tau _3}}𝑑\tau _5\mathrm{exp}[\gamma _2(\tau _5\tau _3)]\sqrt{gϵ^{}}\mathrm{exp}({\displaystyle \frac{\gamma _2\tau _5}{2}})\xi _2^+(\tau _5)`$
$`{\displaystyle _{\tau _3=0}^{\tau _3=\tau }}d\tau _3\mathrm{exp}[\gamma _3(\tau _3\tau )]g{\displaystyle _{\tau _4=0}^{\tau _4=\tau _3}}d\tau _4\mathrm{exp}[\gamma _1(\tau _4\tau _3)]\sqrt{gϵ^{}}\mathrm{exp}({\displaystyle \frac{\gamma _1\tau _4}{2}})\xi _1^+(\tau _4)`$
$`{\displaystyle _{\tau _5=0}^{\tau _5=\tau _3}}d\tau _5\mathrm{exp}[\gamma _2(\tau _5\tau _3)]\sqrt{gϵ^{}}\mathrm{exp}({\displaystyle \frac{\gamma _2\tau _5}{2}})\xi _2^+(\tau _5).`$
However, the average of their product is zero as

$`\alpha _1^{[1/2]}(\tau )\alpha _2^{[1/2]}(\tau )\mathrm{\Delta }\alpha _3^{+[2]}(\tau )`$ (24)
$`=`$ $`g^3ϵϵ^{}{\displaystyle _{\tau _1=0}^{\tau _1=\tau }}{\displaystyle _{\tau _2=0}^{\tau _2=\tau }}{\displaystyle _{\tau _3=0}^{\tau _3=\tau }}{\displaystyle _{\tau _4=0}^{\tau _4=\tau _3}}{\displaystyle _{\tau _5=0}^{\tau _5=\tau _3}}𝑑\tau _1𝑑\tau _2𝑑\tau _3𝑑\tau _4𝑑\tau _5`$ (27)
$`\mathrm{exp}[\gamma _1(\tau _1\tau )]\mathrm{exp}({\displaystyle \frac{\gamma _1\tau _1}{2}})\mathrm{exp}[\gamma _2(\tau _2\tau )]\mathrm{exp}({\displaystyle \frac{\gamma _2\tau _2}{2}})\mathrm{exp}[\gamma _3(\tau _3\tau )]\mathrm{exp}[\gamma _1(\tau _4\tau _3)]\mathrm{exp}({\displaystyle \frac{\gamma _1\tau _4}{2}})`$
$`\mathrm{exp}[\gamma _5(\tau _5\tau _3)]\mathrm{exp}({\displaystyle \frac{\gamma _2\tau _5}{2}})(\xi _1(\tau _1)\xi _2(\tau _2)\xi _1^+(\tau _4)\xi _2^+(\tau _5)\xi _1(\tau _1)\xi _2(\tau _2)\xi _1^+(\tau _4)\xi _2^+(\tau _5))`$
and

$$\xi _1(\tau _1)\xi _2(\tau _2)\xi _1^+(\tau _4)\xi _2^+(\tau _5)=\xi _1(\tau _1)\xi _2(\tau _2)\xi _1^+(\tau _4)\xi _2^+(\tau _5).$$
(28)
Thus, the two noise terms cancel each other and the right hand side of Eq. (28) is zero. Taking such a consideration into account, the moments of the form $`𝒜_1(\tau )𝒜_2(\tau )\mathrm{\Delta }𝒜_3(\tau )`$ are determined by carefully considering the lowest order nonzero terms of their constituent positive-P variables and then finding the average of these variables’ products.

Consider Fig. (2)(a), which contains the lowest order stochastic diagrams for $`\alpha _i(\tau )`$ and $`\alpha _i^+(\tau )`$, where $`i=1,2`$. The first diagram in it represents the initial value terms $`\alpha _i^{(0)}(\tau )`$ and $`\alpha _i^{(0)+}(\tau )`$, which are zero and do not contribute to any moments. The second represents terms containing $`\xi _1,\xi _1^+,\xi _2`$ or $`\xi _2^+`$, which are not necessarily zero and thus may contribute to moments. Fig. 2(b) contains the lowest order stochastic diagrams for $`\alpha _3(\tau )`$ and $`\alpha _3^+(\tau )`$. In it, all terms represented by stochastic diagrams containing initial value lines are zero except for the $`O(g^0)`$ ones. This is so as these terms contain either $`\alpha _i^{(0)}(\tau )\mathrm{or}\alpha _i^{+(0)}(\tau )`$, where $`i=1,2`$, which are both zero. In addition, all $`O(g^0)`$ terms represented by stochastic diagrams in Fig. 2(b) are canceled out by other $`O(g^0)`$ terms. This occurs because $`\mathrm{\Delta }𝒜_3(\tau )`$ appears in the moments considered. Its two components, $`𝒜_3(\tau )`$ and $`𝒜_3(\tau )`$, contain the same $`O(g^0)`$ term and hence their $`O(g^0)`$ terms cancel each other. It follows that the only remaining stochastic diagram in Fig. 2(b), which represents the $`O(g^2)`$ term containing two noise components, denotes the lowest order term in $`\mathrm{\Delta }𝒜_3(\tau )`$ that is not necessarily zero.

The lowest order nonzero terms determined above are now used to calculate $`\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )`$. The lowest order contribution to $`\mathrm{\Delta }\alpha _3(\tau )`$ that is not necessarily zero $`\mathrm{\Delta }\alpha _{3\mathrm{lowest}}(\tau )`$ is

$`\mathrm{\Delta }\alpha _{3\mathrm{lowest}}(\tau )`$ $`=`$ $`{\displaystyle _{\tau _3=0}^{\tau _3=\tau }}𝑑\tau _3\mathrm{exp}[\gamma _3(\tau _3\tau )]g{\displaystyle _{\tau _4=0}^{\tau _4=\tau _3}}𝑑\tau _4\mathrm{exp}[\gamma _1(\tau _4\tau _3)]\sqrt{gϵ}\mathrm{exp}({\displaystyle \frac{\gamma _1\tau _4}{2}})\xi _1(\tau _4)`$ (32)
$`{\displaystyle _{\tau _5=0}^{\tau _5=\tau _3}}𝑑\tau _5\mathrm{exp}[\gamma _2(\tau _5\tau _3)]\sqrt{gϵ}\mathrm{exp}({\displaystyle \frac{\gamma _2\tau _5}{2}})\xi _2(\tau _5)`$
$`{\displaystyle _{\tau _3=0}^{\tau _3=\tau }}d\tau _3\mathrm{exp}[\gamma _3(\tau _3\tau )]g{\displaystyle _{\tau _4=0}^{\tau _4=\tau _3}}d\tau _4\mathrm{exp}[\gamma _1(\tau _4\tau _3)]\sqrt{gϵ}\mathrm{exp}({\displaystyle \frac{\gamma _1\tau _4}{2}})\xi _1(\tau _4)`$
$`{\displaystyle _{\tau _5=0}^{\tau _5=\tau _3}}d\tau _5\mathrm{exp}[\gamma _2(\tau _5\tau _3)]\sqrt{gϵ}\mathrm{exp}({\displaystyle \frac{\gamma _2\tau _5}{2}})\xi _2(\tau _5).`$
When $`\gamma _1=\gamma _2=\gamma _3=\gamma `$, the average of the product of $`\mathrm{\Delta }\alpha _{3\mathrm{lowest}}(\tau )`$ and the lowest order nonzero terms in $`\alpha _1(\tau )`$ and $`\alpha _2(\tau )`$ is approximately equal to $`\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )`$ when $`g1`$ and thus

$`\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )`$ $``$ $`{\displaystyle \frac{ϵ^2g^3\mathrm{exp}[3\gamma \tau ]}{\gamma ^2}}\times \left[{\displaystyle \frac{\mathrm{exp}[\gamma \tau ]}{\gamma }}2\tau {\displaystyle \frac{\mathrm{exp}[\gamma \tau ]}{\gamma }}\right].`$ (33)
As daggered positive-P variables are complex conjugate to undaggered ones on average, $`\alpha _1^+(\tau )\alpha _2^+(\tau )\mathrm{\Delta }\alpha _3^+(\tau )=\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )^{}`$. The other six moments of the form $`𝒜_1(\tau )𝒜_2(\tau )\mathrm{\Delta }𝒜_3(\tau )`$ are zero to $`O(g^3)`$. As does $`\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3^+(\tau )`$, they all have two $`O(g^3)`$ terms that cancel each other. To explain such behaviour in general, the following argument is given. These other six moments can be rewritten as

$$𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau )𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau ),$$
(34)
where it is understood that the moments in which $`𝒜_1(\tau )=\alpha _1(\tau ),𝒜_2(\tau )=\alpha _2(\tau ),𝒜_3(\tau )=\alpha _3(\tau )`$ and $`𝒜_1(\tau )=\alpha _1^+(\tau ),𝒜_2(\tau )=\alpha _2^+(\tau ),𝒜_3(\tau )=\alpha _3^+(\tau )`$ are excluded. All $`O(g^3)`$ terms in the six moments of the form $`𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau )`$ under consideration contain noises in one of the following three forms, $`\xi _1(\tau _a)\xi _2(\tau _b)\xi _1^{}(\tau _c)\xi _2^{}(\tau _d)`$, $`\xi _i(\tau _a)\xi _{3i}^{}(\tau _b)\xi _i(\tau _c)\xi _{3i}(\tau _d)`$ and $`\xi _i^{}(\tau _a)\xi _{3i}(\tau _b)\xi _i^{}(\tau _c)\xi _{3i}^{}(\tau _d)`$, where $`i=1,2`$ and time arguments are dummy variables. All $`O(g^3)`$ terms in the six moments of the form $`𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau )`$ under consideration contain the same noises as their corresponding $`𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau )`$ term. However, in these terms of the form $`𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau )`$ four noise averages from corresponding terms of the form $`𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau )`$ are split into the product of two averages of two noises. For example, $`\alpha _1(\tau )\alpha _2(\tau )\alpha _3^+(\tau )`$ contains noises in the form $`\xi _1(\tau _a)\xi _2(\tau _b)\xi _1^+(\tau _c)\xi _2^+(\tau _d)`$ whilst $`\alpha _1(\tau )\alpha _2(\tau )\alpha _3^+(\tau )`$ contains them in the form $`\xi _1(\tau _a)\xi _2(\tau _b)\xi _1^+(\tau _c)\xi _2^+(\tau _d)`$. Using the formula

$`\xi _1(\tau _a)\xi _2(\tau _b)\xi _3(\tau _c)\xi _4(\tau _d)`$ $`=`$ $`\xi _1(\tau _a)\xi _2(\tau _b)\xi _3(\tau _c)\xi _4(\tau _d)`$ (37)
$`+\xi _1(\tau _a)\xi _3(\tau _c)\xi _2(\tau _b)\xi _4(\tau _d)`$
$`+\xi _1(\tau _a)\xi _4(\tau _d)\xi _2(\tau _b)\xi _3(\tau _c)`$
it can be shown that noise expressions in moments of the form $`𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau )`$ under consideration factorize. In particular, they reduce to the noise expression in the six corresponding term of the form $`𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau )`$. It follows that cancellation occurs between the $`O(g^3)`$ terms in corresponding moments of the form $`𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau )`$ and $`𝒜_1(\tau )𝒜_2(\tau )𝒜_3(\tau )`$ under consideration as the two terms are identical. Consequently, all six moments under consideration are $`O(g^4)`$. They are also typically much smaller than the two $`O(g^3)`$ moments, $`\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )`$ and $`\alpha _1^+(\tau )\alpha _2^+(\tau )\mathrm{\Delta }\alpha _3^+(\tau )`$, as $`g1`$ for realistic systems. To be precise, as all moments are complex quantities, the magnitudes of $`\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )`$ and $`\alpha _1^+(\tau )\alpha _2^+(\tau )\mathrm{\Delta }\alpha _3^+(\tau )`$ are much larger than the magnitudes of the other six moments.

The above results are now more closely related to experiments by considering quadrature phase amplitudes $`X_{i,\theta _i}(\tau )`$. In particular, calculations are performed to determine the in principle experimentally observable third order quadrature phase amplitude moment $`M(\tau )`$, where $`M(\tau )=\mathrm{\Delta }X_{1,\theta _1}(\tau )\mathrm{\Delta }X_{2,\theta _2}(\tau )\mathrm{\Delta }X_{3,\theta _3}(\tau )`$, according to quantum mechanics and SED. In quantum mechanics, quadrature phase amplitudes are expressed in terms of creation and annihilation operators by the equation

$$\widehat{X}_{i,\theta _i}=\frac{\widehat{a}_i\mathrm{exp}(i\theta _i)+\widehat{a}_i^{}\mathrm{exp}(i\theta _i)}{2}.$$
(38)
Using Eq. (38) and operator-positive-P variable correspondences $`\widehat{M}(\tau )_{QM}`$, the value of $`M(\tau )`$ for quantum mechanics can be expressed as

$$\widehat{M}(\tau )_{QM}=\frac{1}{8}\underset{i=1}{\overset{3}{}}\mathrm{\Delta }\alpha _i(\tau )e^{i\theta _i}+\mathrm{\Delta }\alpha _i^+(\tau )e^{i\theta _i}.$$
(39)
Upon expanding the right hand side of Eq. (39), the two lowest order terms in g, $`\mathrm{\Delta }\alpha _1(\tau )\mathrm{\Delta }\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )`$ and $`\mathrm{\Delta }\alpha _1^+(\tau )\mathrm{\Delta }\alpha _2^+(\tau )\mathrm{\Delta }\alpha _3^+(\tau )`$, usually dominate. When they do

$$\widehat{M}(\tau )_{QM}\frac{1}{4}[\mathrm{cos}\mathrm{\Theta }\mathrm{Re}\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )\mathrm{sin}\mathrm{\Theta }\mathrm{Im}\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )],$$
(40)
where $`\mathrm{\Theta }=\theta _1+\theta _2+\theta _3`$. However, when $`\mathrm{cos}\mathrm{\Theta }=0`$ and $`\mathrm{Im}(\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau ))=O(\tau ^4)`$ or when $`\mathrm{sin}\mathrm{\Theta }=0`$ and $`\mathrm{Re}(\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau ))=O(\tau ^4)`$ Eq. (40) is not necessarily true. Such situations can be avoided though because $`\mathrm{\Theta }`$ and $`ϵ`$ are controllable parameters. They are ignored in present considerations. When $`ϵ`$ is real, the $`O(g^3)`$ term in $`\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )`$ is also real and so

$$\widehat{M}(\tau )_{QM}\frac{1}{4}\mathrm{cos}\mathrm{\Theta }\alpha _1(\tau )\alpha _2(\tau )\mathrm{\Delta }\alpha _3(\tau )\frac{ϵ^2g^3\mathrm{cos}\mathrm{\Theta }\mathrm{exp}(3\gamma \tau )}{4\gamma ^2}\left[\frac{\mathrm{exp}(\gamma \tau )}{\gamma }2\tau \frac{\mathrm{exp}(\gamma \tau )}{\gamma }\right].$$
(41)
Thus, Eq. (41) shows that $`\widehat{M}(\tau )_{QM}`$ is cubic in g, within the domain considered, as shown in Fig. 3.

## IV comparison of quantum mechanics and stochastic electrodynamics

This section compares the predictions of quantum mechanics and SED for the intracavity moment $`M(\tau )`$. SED is a semiclassical theory which adds Gaussian white noise to classical electrodynamics. It describes electromagnetic field modes by complex field amplitudes $`\beta `$. For the nondegenerate parametric oscillator, the set of such amplitudes $`\{\beta _1,\beta _2,\beta _3\}`$ evolves via the equations

$`{\displaystyle \frac{\beta _1}{\tau }}`$ $`=`$ $`\gamma _1\beta _1+g\beta _2^{}\beta _3+\sqrt{\gamma _1}\xi _1`$ (42)
$`{\displaystyle \frac{\beta _2}{\tau }}`$ $`=`$ $`\gamma _2\beta _2+g\beta _1^{}\beta _3+\sqrt{\gamma _2}\xi _2`$ (43)
$`{\displaystyle \frac{\beta _3}{\tau }}`$ $`=`$ $`\gamma _3\beta _3g\beta _1\beta _2+\sqrt{\gamma _3}\xi _3,`$ (44)
where the same time variable as in the quantum case is used and the $`\xi ^{}s`$ are independent complex Gaussian white noises with the following correlations

$$\xi _i(\tau _1)\xi _j^{}(\tau _2)=\delta _{ij}\delta (\tau _1\tau _2),$$
(45)
where $`i,j=1,2,3`$. The field amplitudes $`\beta _1,\beta _2\mathrm{and}\beta _3`$ initially have Gaussian fluctuations in their real and imaginary parts of variance 1/4. The only nonzero correlations present in these fluctuations are thus

$$\mathrm{\Delta }\beta _i(0)\mathrm{\Delta }\beta _i^{}(0)=\frac{1}{2},$$
(46)
where $`i=1,2,3`$. Initial conditions are $`\beta _1(0)=\beta _2(0)=0\mathrm{and}\beta _3(0)=ϵ`$.

The SED prediction for the intracavity moment $`M(\tau )`$ is $`M(\tau )_{SED}`$, which is given by the equation

$$M(\tau )_{SED}=\mathrm{\Delta }X_{1,\theta _1}(\tau )\mathrm{\Delta }X_{2,\theta _2}(\tau )\mathrm{\Delta }X_{3\theta _3}(\tau ),$$
(47)
where $`X_{i,\theta _i}(\tau )=\left(\beta _i(\tau )e^{i\theta _i}+\beta _i^{}(\tau )e^{i\theta _i}\right)/2`$. It is calculated using a similar iterative method to the one in Section II, except that noise terms are now treated exactly instead of iteratively. Zeroth order approximations for this iterative method are thus

$`\beta _i^{(0)}(\tau )`$ $`=`$ $`\beta _i(0)\mathrm{exp}(\gamma _i\tau )+{\displaystyle _{\tau _{i=0}}^{\tau _i=\tau }}𝑑\tau _i\mathrm{exp}[\gamma _i(\tau _i\tau )]\sqrt{\gamma _i}\xi _i`$ (48)
$`\beta _3^{(0)}(\tau )`$ $`=`$ $`\beta _3(0)\mathrm{exp}(\gamma _3\tau ),`$ (49)
where $`i=1,2`$. Higher order $`(m+1)^{th}`$ order approximations are

$`\beta _i^{(m+1)}(\tau )`$ $`=`$ $`\beta _i^{(0)}(\tau )+{\displaystyle _{\tau _i=0}^{\tau _i=\tau }}𝑑\tau _i\mathrm{exp}[\gamma _i(\tau _i\tau )]g\beta _{3i}^{(m)}(\tau _i)\beta _3^{(m)}(\tau _i)`$ (50)
$`\beta _3^{(m+1)}(\tau )`$ $`=`$ $`\beta _3^{(0)}(\tau ){\displaystyle _{\tau _3=0}^{\tau _3=\tau }}𝑑\tau _3\mathrm{exp}[\gamma _3(\tau _3\tau )]g\beta _1^{(m)}(\tau _3)\beta _2^{(m)}(\tau _3),`$ (51)
where $`i=1,2`$.

The lowest order nonzero term in g of $`M(\tau )_{SED}`$ is now found using the same method as for the lowest order nonzero term of $`\widehat{M}(\tau )_{QM}`$. Consider the moments of the form $`\mathrm{\Delta }_1(\tau )\mathrm{\Delta }_2(\tau )\mathrm{\Delta }_3(\tau )`$, where $`_n(\tau )`$ is either $`\beta _n(\tau )`$ or $`\beta _n^{}(\tau )`$. The stochastic diagrams required to determine the order of the lowest order nonzero terms of these moments are shown in Fig. 4. Note that noise terms are now $`O(g^0)`$, instead of $`O(g^{\frac{1}{2}})`$ as for quantum mechanics. Using the stochastic diagrams in Fig. 4 it is found that, when $`\gamma _1=\gamma _2=\gamma _3=\gamma `$ and $`g,\tau _f1`$,

$$\mathrm{\Delta }\beta _1(\tau )\mathrm{\Delta }\beta _2(\tau )\mathrm{\Delta }\beta _3^{}(\tau )\frac{g}{12\gamma }[1\mathrm{exp}(3\gamma \tau )]\mathrm{\Delta }\beta _1^{}(\tau )\mathrm{\Delta }\beta _2^{}(\tau )\mathrm{\Delta }\beta _3(\tau )^{}.$$
(52)
The other six moments of the form $`\mathrm{\Delta }_1\mathrm{\Delta }_2\mathrm{\Delta }_3`$ are all $`O(g^2)`$. Thus, $`\mathrm{\Delta }\beta _1(\tau )\mathrm{\Delta }\beta _2(\tau )\mathrm{\Delta }\beta _3^{}(\tau )`$ and $`\mathrm{\Delta }\beta _1^{}(\tau )\mathrm{\Delta }\beta _2^{}(\tau )\mathrm{\Delta }\beta _3(\tau )`$ dominate these other six moments when $`g1`$ and hence

$$M(\tau )_{SED}\frac{1}{4}\mathrm{cos}\mathrm{\Phi }\mathrm{\Delta }\beta _1(\tau )\mathrm{\Delta }\beta _2(\tau )\mathrm{\Delta }\beta _3^{}(\tau )\frac{g\mathrm{cos}\mathrm{\Phi }}{4\gamma }[1\mathrm{exp}(1\gamma \tau )],$$
(53)
where $`\mathrm{\Phi }=\theta _1+\theta _2\theta _3`$, when $`\mathrm{\Phi }0`$. Eq. (53) shows that $`M(\tau )_{SED}`$ is linear in g, as shown in Fig. 3. This is in contrast to the cubic behaviour of $`\widehat{M}(\tau )_{QM}`$. Thus, quantum mechanics and SED predict greatly different values for $`M(\tau )`$ when $`g1`$.

Consideration now is given to the effect of damping strength on the size of the difference between $`M(\tau )_{SED}`$ and $`\widehat{M}(\tau )_{QM}`$. Fig. 5 (a) shows $`\widehat{M}_{QM}`$ and $`M_{SED}`$ as functions of $`\gamma `$ for $`g=0.1,\tau =1`$ and $`ϵ=1`$. It indicates that the difference between them is somewhat sensitive to $`\gamma `$, decreasing exponentially with increasing $`\gamma `$ and quickly approaching zero. However, for the shorter time $`\tau =0.1`$, Fig. 5(b) shows that this difference is not as sensitive to damping. It, approximately, only decreases linearly with increasing $`\gamma `$.

SED and the positive-P representation treat fluctuations very differently, as is evident by comparing noise terms in Eqs (3) and (42). This difference in treatment underlies the differences between the two theories’ results. Firstly, noise terms in the positive-P representation are nonlinear and are scaled by either $`\sqrt{g\alpha _3}`$ or $`\sqrt{g\alpha _3^+}`$, whilst those in SED are linear and are scaled by $`\sqrt{\gamma _i}`$. Secondly, noise terms possess different correlations in the two cases. Thirdly, in quantum mechanics no energy fluctuations occur in the vacuum state, whilst in SED $`\{\beta _i\}`$ fluctuates, as does the total energy. Assuming quantum mechanics is true, in SED fluctuations in the vacuum lead to an overestimate of $`M(\tau )`$ for small g.

## V numerical results

The analytic results for $`M(\tau )_{SED}`$ and $`\widehat{M}(\tau )_{QM}`$ in Sections III and IV only include lowest order nonzero terms. This leaves the sums of all higher order terms as neglected and these may be significant. For this reason, the validity of the analytic approximations are checked by comparison with highly accurate numerical simulation results.

Numerical simulation methods for stochastic differential equations (SDE’s) are both somewhat complex and not widely known. Thus, explanations are given for the numerical technique used to solve the SDE’s in Eqs (3) and (42). Normal ODE techniques such as the Runge-Kutta method cannot be used to solve SDE’s as they contain discontinuous source terms. Instead, a semi-implicit numerical method  is employed. Only its application to Eq. (3) is explained as its application Eq. (42) is similar. Each of the equations in Eq. (3) can be rewritten as

$$\frac{x_i}{\tau }=A_i(𝐱)+\underset{j}{}B_{ij}(𝐱)\zeta _j(\tau ),$$
(54)
where $`x_i`$ is either $`\alpha _i`$ or $`\alpha _i^+`$, for $`i=1,2,3`$, $`𝐱`$ is a vector whose components are $`\{\alpha _i,\alpha _i^+\}`$, $`A_i`$ is the function of $`𝐱`$ formed by the damping and nonlinear terms in the evolution equation for $`x_i`$ and $`b_{ij}`$ is a matrix whose elements are coefficients of the noise terms $`\{\zeta _j\}`$ where $`\zeta _j`$ is either $`\xi _j`$ or $`\xi _j^{}`$, for $`j=1,2`$. The semi-implicit method used determines $`\overline{𝐱}^{(n)}`$, an approximation to $`𝐱`$ at the midpoint of the interval $`(\tau _n,\tau _{n+1})`$. This approximation is found using iteration such that the $`p^{th}`$ order approximation to a component of $`\overline{𝐱}^{(n)}`$ $`\overline{x_i}^{(n)[p]}`$ is given by the equation

$$\overline{x}_i^{(n)[p]}=x_i^{(n)}+\frac{1}{2}\left[\mathrm{\Delta }\tau A_i(\overline{𝐱}^{(n)[p1]})+\underset{j}{}B_{ij}(\overline{𝐱}^{(n)[p1]})\mathrm{\Delta }W_j(\overline{\tau }_n)\right],$$
(55)
where $`x_i^{(n)}`$ is the value of $`x_i`$ at time $`\tau _n`$ , $`\mathrm{\Delta }\tau =\tau _{n1}\tau _n`$, $`\mathrm{\Delta }W_j(\overline{\tau }_n)=\zeta _j^{(n)}(\overline{\tau }_n)\mathrm{\Delta }\tau `$ and $`\overline{\tau }_n`$ is the midpoint of the interval $`(\tau _{n1},\tau _n)`$. The zeroth order approximation to $`\overline{x}_i^{(n)}`$ is given by the equation $`\overline{x}_i^{(n)[0]}=x_i^{(n)}`$. The approximation to $`\overline{𝐱}^{(n)}`$ calculated is then used to generate $`\mathrm{\Delta }x_i^{(n)}`$, an approximation to the change in $`x_i`$ over the interval $`(\tau _n,\tau _{n+1})`$. This is done by solving the equation

$$\mathrm{\Delta }x_i^{(n)}=A_i(\overline{𝐱}^{(n)})\mathrm{\Delta }\tau _n+\underset{j}{}B_{ij}(\overline{𝐱}^{(n)})\mathrm{\Delta }W_j(\overline{\tau }_n).$$
(56)
Repeated use of Eq. (56) determines $`x_i^{(n)}`$ for successively later and later times and thus solves Eq. (54). Two of the most important parameters used in the numerical simulations are the step size and the number of stochastic paths that are averaged over. The former is always 0.0025 and the latter is $`O(10^6)`$ for most simulations. However, large sampling errors necessitated averaging over $`O(10^7)`$ paths for $`g=0.1`$ SED simulations.

Results from the numerical and analytic simulations of $`M(\tau )_{SED}`$ and $`\widehat{M}(\tau )_{QM}`$ over a range of $`g`$ and $`N`$ values, where $`N`$ is the average initial number of pump photons ($`N=|ϵ^2|`$), are shown in Figs 6-8. In all cases $`\theta _1=\theta _2=\theta _3=0`$ and relative numerical errors are small. All $`g=0.1`$ analytic results are in agreement with their numerical counterparts. However, $`g=1`$ analytic results for N=1 and N=10 are not. This disagreement is explained by noting that the analytic results are only necessarily valid when $`g1`$.

A number of qualitative trends can be seen in Figs 6-8. In Figs 6 and 7 (N=1 and N=10) the results of SED and quantum mechanics are so distinct that they have different signs, with those of quantum mechanics being negative and those of SED being positive. This trend only holds for short times ($`\tau <0.07`$) in Figs 8(a) and (b) (N=100). For longer times, SED and quantum mechanics predict the same sign. This trait is consistent with the fact that Figs 8(a) and (b) show results for the largest number of photons in the pump mode. SED and quantum mechanics are at their most classical level for this case and thus might be expected to differ the least. For constant N Figs 6-8 also show that as $`g`$ is decreased the results of quantum mechanics and SED become more similar. This occurs because lower $`g`$ values are associated with larger damping to nonlinear coupling ratios and therefore move SED and quantum mechanics closer to the classical domain.

## VI External moments

Thus far, only intracavity fields have been considered. However, it is the external fields that leak out of a cavity that are observed. In realistic systems, intracavity photons are transmitted through imperfect mirrors into the external environment where they are detected. Thus, an external field analogue of $`M(\tau )`$, $`M^{(E)}(\tau _s,\tau _f)`$, where $`\tau _s`$ and $`\tau _f`$ are initial and final measurement times, is calculated according to quantum mechanics and SED to consider what is actually observed in the laboratory.

The first step in calculating the external moment $`\widehat{M}^{(E)}(\tau )_{QM}`$ for quantum mechanics is defining the external quadrature phase amplitudes constituting it. This is done within the the context of homodyne detection as quadrature phase amplitudes are commonly measured using it. A schematic diagram for balanced homodyne detection is shown in Fig. 9. An external signal field flux $`\widehat{\mathrm{\Phi }}_{iOUT}`$, where $`i=1,2,3`$, and a local oscillator field flux $`E_i`$ are incident on a 50-50 beam splitter BS. An external local oscillator phase variable is represented by $`\overline{\theta }_i`$. The two field fluxes combine and are detected by two photodiodes $`D_{+i}`$ and $`D_i`$. The detected photocurrents are then converted to amplified electrical currents whose difference is found. An external quadrature phase amplitude for quantum mechanics $`\widehat{X}_{i,\overline{\theta }_i}^{(E)}`$ is defined as this difference yielding, when $`E_i`$ is real,

$$\widehat{X}_{i,\overline{\theta }_i}^{(E)}(\tau )=\frac{eA\eta _iE_i(\widehat{\mathrm{\Phi }}_{iOUT}(\tau )e^{i\overline{\theta }_i}+\widehat{\mathrm{\Phi }}_{iOUT}^{}{}_{}{}^{}(\tau )e^{i\overline{\theta }_i})}{2},$$
(57)
where e is the magnitude of the charge of an electron, A is an amplification factor and $`\eta _i`$ is a detector efficiency factor for both detectors associated with external field modes denoted by $`i`$. In realistic experiments detection occurs over a finite period of time and thus

$$_{\tau =\tau _s}^{\tau =\tau _f}\frac{d\tau }{\mathrm{\Gamma }}\widehat{X}_{i,\overline{\theta }_i}^{(E)}(\tau )$$
(58)
corresponds to what is observed. Only the $`\overline{\theta }_i=0`$ case is considered. Thus, an external moment analogue of the intracavity moment $`\widehat{M}(\tau )_{QM}`$ can be defined as

$`\widehat{M}^{(E)}(\tau _s,\tau _f)_{QM}`$ $`=`$ $`\mathrm{\Gamma }^3{\displaystyle \underset{i=1}{\overset{3}{}}}{\displaystyle _{\tau _i=\tau _s}^{\tau _i=\tau _f}}𝑑\tau _i\mathrm{\Delta }\widehat{X}_{i,\overline{\theta }_i}^{(E)}(\tau _i)`$ (59)
$`=`$ $`\mathrm{\Gamma }^3{\displaystyle _{\tau _1=\tau _s}^{\tau _1=\tau _f}}{\displaystyle _{\tau _2=\tau _s}^{\tau _2=\tau _f}}{\displaystyle _{\tau _3=\tau _s}^{\tau _3=\tau _f}}{\displaystyle \underset{i=1}{\overset{3}{}}}\mathrm{\Delta }\widehat{X}_{i,\overline{\theta }_i}^{(E)}(\tau _i).`$ (60)
To calculate $`\widehat{M}^{(E)}(\tau _s,\tau _f)_{QM}`$, the relation between the unknown external output fields that define it and known intracavity fields needs to be ascertained. Gardiner and Collett  have formulated an input-output theory which relates the two via the equation

$$\widehat{\mathrm{\Phi }}_{iOUT}(\tau )=\sqrt{2\mathrm{\Gamma }}\widehat{a}_i(\tau )+\widehat{\mathrm{\Phi }}_{iIN}(\tau ),$$
(61)
where $`\widehat{\mathrm{\Phi }}_{iIN}(\tau )`$ is the input field flux associated with intracavity mode $`i`$. All input fields are assumed to be in vacuum states. This allows the use of Eq. (5.3) from , which can be expressed as, in this paper’s notation,

$`\widehat{\mathrm{\Phi }}_{iOUT}^{}(\tau _1)\widehat{\mathrm{\Phi }}_{iOUT}^{}(\tau _2)\mathrm{}\widehat{\mathrm{\Phi }}_{iOUT}^{}(\tau _n)\widehat{\mathrm{\Phi }}_{iOUT}(\tau _{n+1}^{})\mathrm{}\widehat{\mathrm{\Phi }}_{iOUT}(\tau _m^{})`$ (62)
$`=(2\mathrm{\Gamma })^{m/2}\stackrel{~}{T}[\widehat{a}_i^{}(\tau _1)\widehat{a}_i^{}(\tau _2)\mathrm{}\widehat{a}_i^{}(\tau _n)]T[\widehat{a}_i(\tau _{n+1}^{})\mathrm{}\widehat{a}_i(\tau _m^{})]`$ , (63)
where $`\stackrel{~}{T}`$ and $`T`$ are time anti-ordering and time ordering operators respectively. Using Eq. (57), the integrand of Eq. (59) can be expressed in terms of $`\widehat{\mathrm{\Phi }}_{iOUT}`$ and $`\widehat{\mathrm{\Phi }}_{iOUT}^{}`$. It can then be expressed in terms of particular $`\widehat{a}_i(\tau _i)`$ and $`\widehat{a}_{i}^{}{}_{}{}^{}(\tau _i)`$ averages using Eq. (62). In turn, these averages are equivalent to the eight positive-P averages of the form $`\mathrm{\Delta }𝒜_1(\tau _1)\mathrm{\Delta }𝒜_2(\tau _2)\mathrm{\Delta }𝒜_3(\tau _3)`$, where $`𝒜_i`$ is either $`\alpha _i`$ or $`\alpha _i^+`$. As was determined in Section III, two of these averages, $`\mathrm{\Delta }\alpha _1(\tau _1)\mathrm{\Delta }\alpha _2(\tau _2)\mathrm{\Delta }\alpha _3(\tau _3)`$ and $`\mathrm{\Delta }\alpha _1^+(\tau _1)\mathrm{\Delta }\alpha _2^+(\tau _2)\mathrm{\Delta }\alpha _3^+(\tau _3)`$, are of lower order in $`g`$ than the others and hence dominate when $`g1`$. Thus, the external field moment of quantum mechanics $`\widehat{M}^{(E)}(\tau _s,\tau _f)_{QM}`$ can be expressed as, when $`\overline{\theta }_i=0`$, where $`i=1,2,3`$,

$`\widehat{M}^{(E)}(\tau _s,\tau _f)_{QM}`$ $``$ $`{\displaystyle \frac{\sqrt{2}(eA\eta E)^3\mathrm{\Gamma }^{3/2}}{4}}{\displaystyle _{\tau _1=\tau _s}^{\tau _1=\tau _f}}{\displaystyle _{\tau _2=\tau _s}^{\tau _2=\tau _f}}{\displaystyle _{\tau _3=\tau _s}^{\tau _3=\tau _f}}𝑑\tau _1𝑑\tau _2𝑑\tau _3`$ (65)
$`\mathrm{\Delta }\alpha _1(\tau _1)\mathrm{\Delta }\alpha _2(\tau _2)\mathrm{\Delta }\alpha _3(\tau _3)+\mathrm{\Delta }\alpha _1^+(\tau _1)\mathrm{\Delta }\alpha _2^+(\tau _2)\mathrm{\Delta }\alpha _3^+(\tau _3),`$
where $`\eta =\eta _i`$ and $`E=E_i`$, where $`i=1,2,3`$. To simplify the algebra only the $`\tau _s=0`$ case is investigated, so that only the moment $`\widehat{M}^{(E)}(\tau _f)_{QM}(\widehat{M}^{(E)}(0,\tau _f)_{QM})`$ is considered. External fields are only considered for small times ($`\tau _f<<1`$), and so, to a given order in $`g`$, $`\widehat{M}^{(E)}(\tau _f)_{QM}`$’s lowest nonzero order term in $`\tau _f`$ dominates. Hence $`\widehat{M}^{(E)}(\tau _f)_{QM}`$ can be approximated by its lowest nonzero order term in both $`g`$ and $`\tau _f`$. Thus,

$$\widehat{M}^{(E)}(\tau _f)_{QM}\frac{\sqrt{2}}{48}g^3ϵ^2(eA\eta E)^3\mathrm{\Gamma }^{3/2}\tau _f^6.$$
(66)
The SED external moment $`M^{(E)}(\tau _f)_{SED}`$ is now calculated. It is given by the same expression as $`\widehat{M}(\tau _f)_{QM}`$, the right hand side of Eq. (59) (when $`\tau _s=0`$), except that the quadrature phase amplitude operator $`\widehat{X}_{i\overline{\theta }_i}^{(E)}(\tau )`$, is replaced by its SED c-number analogue. This external SED c-number quadrature phase amplitude is defined as, when $`E_i`$ is real,

$$X_{i,\overline{\theta }_i}^{(E)}(\tau )=\frac{eA\eta _iE_i(\beta _{iOUT}(\tau )e^{i\overline{\theta }_i}+\beta _{iOUT}^{}(\tau )e^{i\overline{\theta }_i})}{2},$$
(67)
where $`\beta _{iOUT}(\tau )`$ is the output field flux associated with the intracavity field denoted by $`i`$. In analogy with Eq. (61), it is assumed that the SED input-output relation is

$$\beta _{iOUT}(\tau )=\sqrt{2\mathrm{\Gamma }}\beta _i(\tau )+\beta _{iIN}(\tau ),$$
(68)
where $`\beta _{iIN}(\tau )`$ is the input field flux for the intracavity mode $`i`$. When all input fields are in vacuum states, as is the case, $`\beta _{iIN}(\tau )`$ is a Gaussian white noise with a self correlation characterized by

$$\beta _{iIN}(\tau _i)\beta _{iIN}^{}(\tau _i^{})=\frac{\delta (\tau _i\tau _i^{})}{2\mathrm{\Gamma }}.$$
(69)
A calculation analogous to the quantum mechanical one earlier in this section can be performed using Eqs (67) and (68) to obtain an expression for $`M^{(E)}(\tau _f)_{SED}`$ in terms of particular intracavity averages. When lowest order nonzero approximations to these averages are considered, the following result is obtained when $`\overline{\theta }_i`$ for $`i=1,2,3`$, and $`g,\tau _f1`$,

$$M^{(E)}(\tau _f)_{SED}\frac{\sqrt{2}}{16}g\tau _f^4(eA\eta E)^3\mathrm{\Gamma }^{3/2}.$$
(70)
Upon comparing Eq. (70) to result of quantum mechanics in Eq. (66), it is seen that the leading order term in $`g`$ in Eq. (66) is $`O(g^3)`$ whilst in Eq. (70) it is $`O(g)`$. Hence, as was the case for the intracavity moment, quantum mechanics and SED predict significantly different results for the observable external field moment $`M^{(E)}(\tau _f)`$.

## VII Signal to noise ratio

In actual experiments, only finite samples of results are obtained, as opposed to infinite ones. Hence, in practice the population means considered thus far are estimated from sample means. These sample means fluctuate from sample to sample and thus have signal to noise ratios, which are now determined for small times ($`\tau _f<<1`$). This paper focuses on differences between quantum mechanics and SED. Thus, a calculation is performed of the signal to noise ratio of the difference between the two theories’ external sample moments. First, the noise of the external sample moment in quantum mechanics is determined. It is then assumed that the noise of the external sample moment of SED is the same. Noise results are combined with the external moment results of Section VI to produce $`S(\tau _f)`$, the signal to noise ratio of the difference between the two theories’ external sample moments. This quantity $`S(\tau _f)`$ is given by

$`S(\tau _f)`$ $`=`$ $`{\displaystyle \frac{|m^{(E)}(\tau _f)_{SED}\widehat{m}^{(E)}(\tau _f)_{QM}|}{\sqrt{s^2(m^{(E)}(\tau _f)_{SED})+s^2(\widehat{m}^{(E)}(\tau _f)_{QM})}}}`$ (71)
$`=`$ $`{\displaystyle \frac{\sqrt{n1}|M^{(E)}(\tau _f)_{SED}\widehat{M}^{(E)}(\tau _f)_{QM}|}{\sqrt{2}\sigma (\widehat{M}^{(E)}(\tau _f))}},`$ (72)
where $`n`$ is the number of observations in the sample considered, $`m^{(E)}(\tau _f)_{SED}`$ and $`\widehat{m}^{(E)}(\tau _f)_{QM}`$ are sample averages of $`M^{(E)}(\tau _f)`$ according to SED and quantum mechanics respectively, and $`s^2(A)`$ denotes the sample variance of A. The only significant unknown quantity on the right hand side of Eq. (71) is $`\sigma (\widehat{M}^{(E)}(\tau _f))`$, which is now determined. Expressing $`\sigma (\widehat{M}^{(E)})`$ explicitly yields

$$\sigma (\widehat{M}^{(E)}(\tau _f))=\sqrt{\widehat{M}^{(E)}(\tau _f)^2_{QM}\widehat{M}^{(E)}(\tau _f)_{QM}^2}.$$
(73)
The moment $`\widehat{M}^{(E)}(\tau _f)_{QM}`$ was determined in Section VI and so $`\widehat{M}^{(E)}(\tau _f)^2_{QM}`$ is now calculated. In the calculation that follows only the $`\overline{\theta }_i=0`$, where $`i=1,2,3`$, and $`g,\tau _f1`$ cases are considered.

The moment $`\widehat{M}^{(E)}(\tau _f)^2_{QM}`$ can be expressed in terms of external quadrature phase operators as

$`\widehat{M}^{(E)}(\tau _f)^2_{QM}`$ $`=`$ $`[\mathrm{\Gamma }^3{\displaystyle \underset{i=1}{\overset{3}{}}}{\displaystyle _{\tau _i=0}^{\tau _i=\tau _f}}𝑑\tau _i\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i)]^2`$ (74)
$`=`$ $`\mathrm{\Gamma }^6{\displaystyle _{\tau _1=0}^{\tau _1=\tau _f}}{\displaystyle _{\tau _1^{}=0}^{\tau _1^{}=\tau _f}}{\displaystyle _{\tau _2=0}^{\tau _2=\tau _f}}{\displaystyle _{\tau _2^{}=0}^{\tau _2^{}=\tau _f}}{\displaystyle _{\tau _3=0}^{\tau _3=\tau _f}}{\displaystyle _{\tau _3^{}=0}^{\tau _3^{}=\tau _f}}𝑑\tau _1𝑑\tau _1^{}𝑑\tau _2𝑑\tau _2^{}𝑑\tau _3𝑑\tau _3^{}{\displaystyle \underset{i=1}{\overset{3}{}}}\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i)\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i^{}),`$ (75)
where $`\widehat{X}_i^{(E)}(\tau _i)=\widehat{X}_{i\overline{\theta }_i=0}^{(E)}(\tau _i).`$ The integrand of Eq. (74), which is denoted by K, can be expressed as

$$K=\underset{i=1}{\overset{3}{}}\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i)\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i^{})+f(\tau _i,\tau _i^{}),$$
(76)
where $`f(\tau _i,\tau _i^{})`$ is a function which includes terms resulting from coupling between modes. These coupling terms vanish when $`g=0`$ and thus are at least O(g). It follows that $`K`$ can be be re-expressed as

$$K=\underset{i=1}{\overset{3}{}}\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i)\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i^{})+O(g).$$
(77)
The moment $`\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i)\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i^{})`$ is now calculated using a normally ordered approach that has been previously employed to solve similar problems ,. This method expresses $`\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i)\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i^{})`$ in terms of normally ordered photocurrent averages and then determines these averages. It first defines $`\widehat{X}_i^{(E)}(\tau _a)`$, where $`\tau _a`$ is any $`\tau `$ variable, as the difference between the amplified electrical currents, $`\widehat{X}_{+i}^{(E)}(\tau _a)`$ and $`\widehat{X}_i^{(E)}(\tau _a)`$, produced by the photocurrents detected at the detectors $`D_{+i}`$ and $`D_i`$ in Fig. 9 in Section VI. Using this definition ($`\widehat{X}_i^{(E)}(\tau _a)=\widehat{X}_{+i}^{(E)}(\tau _a)\widehat{X}_i^{(E)}(\tau _a)`$), $`\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i)\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i^{})`$ can be expressed as

$$\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i)\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i^{})=(\mathrm{\Delta }\widehat{X}_{+i}^{(E)}(\tau _i)\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i))(\mathrm{\Delta }\widehat{X}_{+i}^{(E)}(\tau _i^{})\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i^{})).$$
(78)
Upon expansion, the right hand side of Eq. (78) contains two types of terms, those of the form $`\widehat{X}_{Ci}^{(E)}(\tau _a)`$, where C is either + or -, and those of the form $`\widehat{X}_{Ci}^{(E)}(\tau _i)\widehat{X}_{Di}^{(E)}(\tau _i^{})`$, where D is either + or -. Terms of the form $`\widehat{X}_{Ci}^{(E)}(\tau _a)`$ are given by the equation

$$\widehat{X}_{Ci}^{(E)}(\tau _a)=_{\mathrm{}}^+\mathrm{}\frac{ds_1}{\mathrm{\Gamma }}G_{Ci}^{(1)}(s_1)J^0(\tau _as_1),$$
(79)
where $`G_{Ci}^{(1)}(s_1)`$ is a first order Glauber correlation function and $`J^0(\tau _as_1)`$ is an electrical current pulse produced by a single photodetection event. In following previous work , square electrical current pulses of the form

$`J^0\left(ab\right)=\{\begin{array}{cc}Ae\mathrm{\Gamma }/\tau _d\hfill & bab+\tau _d\hfill \\ 0\hfill & a<b\mathrm{and}a>b+\tau _d\hfill \end{array}`$ (82)
are considered in the limit of $`\tau _d0`$, which is taken at some appropriate later stage of the calculation. The Glauber correlation function $`G_{Ci}^{(1)}(s_1)`$ can be expressed as a power series in $`g`$ and $`s_1`$ and thus as $`_{m,n=0}^\mathrm{},\mathrm{}c_{mn}g^ms_1^n`$. Due to the form of $`J^0(\tau _as_1)`$, when $`\tau _a<<1`$, as is being assumed, only photodetection events at small times $`s_1`$, contribute to $`\widehat{X}_{Ci}^{(E)}(\tau _a)`$. This fact, coupled with the knowledge that only the $`g<<1`$ case is considered, means that the $`n=m=0`$ term in the power series for $`G_{Ci}^{(1)}(s_1)`$ dominates when $`J^0(\tau _as_1)`$ is nonzero. Hence, upon calculating this dominant term by expressing $`\widehat{\mathrm{\Phi }}_i`$ and $`\widehat{\mathrm{\Phi }}_i^{}`$ in terms of intracavity field operators, in the limit of large local oscillator amplitude,

$$G_{Ci}^{(1)}(s_1)\frac{\eta _{Ci}}{2}E_i^2,$$
(83)
where $`\eta _{Ci}`$ is a detector efficiency factor for the photodetector $`D_{Ci}`$. It follows that

$$\widehat{X}_{Ci}^{(E)}(\tau _a)\frac{\eta _{Ci}E_i^2Ae}{2}.$$
(84)
Terms of the form $`\widehat{X}_{Ci}^{(E)}(\tau _i)\widehat{X}_{Di}^{(E)}(\tau _i^{})`$ in Eq. (78) can be expressed as

$`\widehat{X}_{Ci}^{(E)}(\tau _i)\widehat{X}_{Di}^{(E)}(\tau _i^{})`$ $`=`$ $`\delta _{CD}{\displaystyle _{\mathrm{}}^+\mathrm{}}{\displaystyle \frac{ds_1}{\mathrm{\Gamma }}}G_{Ci}^{(1)}(s_1)J^{(0)}(\tau _is_1)J^{(0)}(\tau _i^{}s_1)`$ (86)
$`+{\displaystyle _{\mathrm{}}^+\mathrm{}}{\displaystyle _{\mathrm{}}^+\mathrm{}}{\displaystyle \frac{ds_1ds_2}{\mathrm{\Gamma }^2}}G_{C,Di}^{(2)}(s_1,s_2)J^{(0)}(\tau _is_1)J^{(0)}(\tau _i^{}s_2),`$
where $`G_{C,Di}^{(2)}(s_1,s_2)`$ is a second order Glauber correlation function and $`\delta _{C,D}`$ is one when C and D are the same and zero otherwise. In the limit of $`\tau _d0`$,

$$_{\mathrm{}}^+\mathrm{}\frac{ds_1}{\mathrm{\Gamma }}G_{Ci}^{(1)}(s_1)J^{(0)}(\tau _is_1)J^{(0)}(\tau _i^{}s_1)\frac{(eAE_i)^2\eta _{Ci}\delta (\tau _i\tau _i^{})\mathrm{\Gamma }}{2}$$
(87)
to leading nonzero order in $`g,\tau _i`$ and $`\tau _i^{}`$. It is of equal order in $`g`$ and lower order in $`\tau _i`$ and $`\tau _i^{}`$ than the second term in Eq. (86) and hence is much larger than this second term when it is nonzero as the $`\tau _i,\tau _i^{}<<1`$ case is being considered. Thus

$$\widehat{X}_{Ci}^{(E)}(\tau _i)\widehat{X}_{Di}^{(E)}(\tau _i^{})\frac{\delta _{CD}(AeE_i)^2\eta _{Ci}\delta (\tau _i\tau _i^{})\mathrm{\Gamma }}{2}.$$
(88)
From Eqs (84) and (88) it can be seen that the single integral terms in $`\widehat{X}_{+i}^{(E)}(\tau _i)\widehat{X}_{+i}^{(E)}(\tau _i^{})`$ and $`\widehat{X}_i^{(E)}(\tau _i)\widehat{X}_i^{(E)}(\tau _i^{})`$ are of the same order in $`g`$ and lower order in $`\tau _i`$ and $`\tau _i^{}`$ than any other terms contributing to $`\widehat{X}_i^{(E)}(\tau _i)\widehat{X}_i^{(E)}(\tau _i^{})`$ and hence dominate. It follows that

$$\widehat{X}_i^{(E)}(\tau _i)\widehat{X}_i^{(E)}(\tau _i^{})(AeE)^2\eta _i\delta (\tau _i\tau _i^{})\mathrm{\Gamma },$$
(89)
where $`\eta _i=\eta _{Ci}=\eta _{Di}`$ and $`E=E_i`$, where $`i=1,2,3`$. As right hand side of Eq. (89) is $`O(g^0)`$, $`_{i=1}^3\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i)\mathrm{\Delta }\widehat{X}_i^{(E)}(\tau _i^{})`$ is also $`O(g^0)`$ and hence from Eq. (77),

$$K((AeE)^2\mathrm{\Gamma })^3\underset{i=1}{\overset{3}{}}\eta _i\delta (\tau _i\tau _i^{}).$$
(90)
Substituting this approximation for K into Eq. (74) yields

$`\widehat{M}^{(E)}(\tau _f)^2`$ $``$ $`\mathrm{\Gamma }^6((eAE)^2\eta \mathrm{\Gamma })^3{\displaystyle \underset{i=1}{\overset{3}{}}}{\displaystyle _{\tau _i=0}^{\tau _i=\tau _f}}{\displaystyle _{\tau _i^{}}^{\tau _i^{}=\tau _f}}𝑑\tau _i𝑑\tau _i^{}\delta (\tau _i\tau _i^{})`$ (91)
$`=`$ $`\left[{\displaystyle \frac{(eAE)^2\eta \tau _f}{\mathrm{\Gamma }}}\right]^3,`$ (92)
where $`\eta =\eta _i`$, where $`i=1,2,3`$. Thus

$$\sigma (\widehat{M}^{(E)}(\tau _f))\left[\frac{(eAE)^2\eta \tau _f}{\mathrm{\Gamma }}\right]^{3/2}.$$
(93)
Hence, the signal to noise ratio of the difference between the external sample moments of quantum mechanics and SED is

$$S(\tau _f)\frac{\sqrt{n1}\eta ^{3/2}g\tau _f^{5/2}}{16}.$$
(94)
## VIII Realistic systems

Realistic parameter values are now considered to determine if the theoretical difference between SED and quantum mechanics could be observed experimentally. In particular, the signal to noise ratio of the difference between the sample moments of quantum mechanics and SED $`S`$ is calculated using realistic parameter values for nondegenerate parametric oscillators containing the commonly used crystals, silver gallium selinide ($`\mathrm{AgGaSe}_2`$) and potassium titanyl phosphate (KTP). The non-linear interaction strength G for parametric down conversion is given by

$$Gd_{eff}\sqrt{\frac{2\mathrm{}\omega _1\omega _2\omega _3}{ϵ_0V}}\frac{l}{L},$$
(95)
where V is the cavity volume, l is the crystal length and L the cavity length. Cavity and crystal length values of 10cm are chosen. The cavity volume V is given by the formula $`V=\pi \mathrm{\Omega }^2L`$, where $`\mathrm{\Omega }`$ is the spot size. This volume is minimized in order to maximize G and thus the external difference between quantum mechanics and SED. It is assumed that the damping constant used to scale time $`\mathrm{\Gamma }`$ equals the unscaled damping constant for each mode $`\mathrm{\Gamma }_i`$ ($`\mathrm{\Gamma }=\mathrm{\Gamma }_i`$). This common damping constant $`\mathrm{\Gamma }`$ is calculated from the formula $`\mathrm{\Gamma }=T\times \frac{c}{2L}`$, where c is the speed of light and T is a mirror transmission coefficient. A $`T`$ value of $`T=0.01`$ is used. Using the above information, Table I shows realistic parameter values for $`d_{eff}`$, $`V`$, $`G`$, pump, signal and idler wavelengths, and resulting $`g`$ and $`\mathrm{\Gamma }`$ values. Results for $`\widehat{M}^{(E)}_{QM}`$ and $`M^{(E)}_{SED}`$ are obtained using Eqs (66) and (70) for when $`\eta =1`$, $`\tau _f=0.1`$, $`E=10^9s^{1/2}`$, $`A=1/e`$ and $`ϵ=10^3`$. These are displayed in Table II, which shows that the external results of quantum mechanics and SED differ greatly. Due to local oscillator amplification, they are also macroscopically distinct with respect to photon number, even though the initial number of intracavity photons is small on average. Another appealing feature of the difference between the two theories is that detector efficiencies approach one as photodiodes as opposed to photomultipliers are used for detection. Thus, no fair sampling assumptions need to be made.

The question remains of whether or not the population difference between SED and quantum mechanics could be reliably observed in a finite sample of results. To answer it, $`S`$ is now considered. Figs 10 and 11 show graphs of $`S`$ versus sample size $`n`$ for $`\mathrm{AgGaSe}_2`$ and KTP for the same parameter values as used in the last paragraph. These show reasonable $`S`$ values and indicate that large sample sizes must be obtained to produce a signal to noise ratio of one, the smallest signal to noise ratio required to clearly observe the signal. In particular, sample sizes of $`1.8\times 10^{13}`$ (KTP) and $`3.7\times 10^{11}`$ ($`\mathrm{AgGaSe}_2`$) need to be obtained to generate a signal to noise ratio of one. An individual observation takes a time of the order $`t=\tau _f/\mathrm{\Gamma }=6.7\times 10^9s`$ and so, assuming minimal time delay between measurements, $`1.8\times 10^{13}`$ observations would take about 33 hours and $`3.7\times 10^{11}`$ observations about 41 minutes. It is conceivable that measurements could be taken over both times. Furthermore, as the signal to noise ratio scales as $`\frac{1}{g}`$, higher $`g`$ materials would enable the difference to be observed even more readily.

## IX discussion

It has been shown that there exist a significant, potentially experimentally observable, difference between quantum mechanics and SED. Due to local oscillator amplification, this difference can involve macroscopically distinct external fields for the two theories. Thus, it can be considered macroscopic if it is legitimate to include the local oscillators as part of the system and not as external measuring apparatuses. The difference is also potentially experimentally observable, as a realistic system and state are considered and is present at realistic parameter values. The system is practical as parametric oscillators and balanced homodyne detection are widely used, and damping is included. The state is realistic as the initial intracavity coherent state can be approximated well by a laser. It follows that the difference can be seen as providing the basis for an experimentally achievable macroscopic test of quantum mechanics against one local hidden variable theory (SED). Such a test is significant as all experimental tests of quantum mechanics against local hidden variable theories to date have been microscopic. It is true that many macroscopic tests have been proposed, but most of them consider highly idealized states or systems that are not currently able to be experimentally implemented. In particular, many of them do not consider damping, even though it is known to rapidly destroy the correlations of quantum mechanics present in Schroedinger cat  and other entangled states. The calculations in this paper do include damping and show that the difference between SED and quantum mechanics is not overly sensitive to it. Most importantly, it remains for realistic damping values. The test proposed in this paper can be seen as being in the novel and largely unexplored domain of macroscopic experimental tests of quantum mechanics.

Even if the local oscillators are not included as part of the system investigated, the external difference between quantum mechanics and SED is still at least mesoscopic as average initial pump photon numbers up to $`10^6`$ are considered. From this perspective, the difference is still distinct from many earlier microscopic ones known to exist between quantum mechanics and all local hidden variable theories. It is also, perhaps, more surprising than some of them as it occurs in a larger particle number system.

Two noteworthy features of the external difference between quantum mechanics and SED are that it involves continuous variables and high efficiency detection. That it involves continuous variables is significant because most previous differences between quantum mechanics and local hidden variable theories have involved discrete ones. Furthermore, it is, perhaps, more surprising that a difference between quantum mechanics and a local hidden variable theory can be found for continuous variables as continuous variables are more closely related to classical ones (which are all continuous) than discrete ones. Low detector efficiency forms the basis of a significant loophole in most tests between quantum mechanics and SED to date . The use of photodiodes for detection in the scheme discussed means that such a loophole is avoided.

The calculations in Section VIII show it is difficult to observe the external difference between quantum mechanics and SED. This is mainly a result of small experimental nonlinearities. They cause few signal and idler photons to be created and thus the experimental signal is weak relative to its noise. For small enough measurement samples, SED results cannot be clearly distinguished from those of quantum mechanics. This fact is consistent with the knowledge that SED reproduces many features of quantum mechanics. However, it is a distinct theory and does differ from quantum mechanics in particular cases, as this paper has shown.

The external difference between quantum mechanics and SED would be easier to observe if larger nonlinear coupling constants were used. These could be achieved by using organic nonlinear crystals such as N-(4-nitrophenyl)-L-prolinol (NPP) . However, phase matching would be difficult with such crystals. In addition, they are typically only transparent within a small frequency range. Alternatively, higher nonlinearities could be achieved by using Josephson-parametric amplifiers , which can have even larger nonlinearities than organic nonlinear crystals. Another possibility, in the area of atom optics, is to to utilize BEC nonlinear effects, in which atom-molecule coupling is induced through photon-associaton .

To conclude, this paper compared particular moments of quantum mechanics to those of SED for the nondegenerate parametric oscillator. Both internal and external moments were considered and an analytic iterative technique showed them both to be cubic in the system’s nonlinear coupling constant for quantum mechanics and linear for SED. Numerical simulations were performed to check the approximate intracavity analytic result and were in agreement with them when the system’s nonlinear coupling constant was much less than one. Realistic parameter values were considered and it was shown that the external sample difference between SED and quantum mechanics had a small signal to noise ratio in typical parametric oscillators. The presence of intense local oscillators means that the results could be seen as providing the basis for a macroscopic experimental test of quantum mechanics against SED.

## acknowledgements

DTP would like to thank Professor Gerard Milburn, Dr Howard Wiseman, Dr Karen Kheruntsyan, Professor Brian Orr, Michael Gagen and Cynthia Freeman for their assistance with the paper. He would also like to thank The University of Queensland for its financial support. WJM acknowledges support from the Australian Research Council.