# Lasing at half the Josephson frequency with exponentially long coherence times

###### Abstract

We describe a superconducting device capable of producing laser light in the visible range at half the Josephson generation frequency, with the optical phase of the light locked to the superconducting phase difference. An earlier proposed device, the so called “half-Josephson laser” [Phys. Rev. Lett. 107, 073901 (2011)], cannot provide long coherence times, because of spontaneous switchings between the emitter states. To circumvent this we consider emitters driving an optical resonator mode. We derive a general model that captures essential physics of such devices while not depending on specific microscopic details. We find the conditions under which the coherence times are exponentially long, thus surpassing the fundamental limitation on the coherence times of common lasers. For this we study the noise in the device. In particular, we are interested in the rate of large fluctuations of the light field in the limit where the typical fluctuations are small. The large fluctuations are responsible for switching of the laser between stable states of radiation and therefore determine the coherence time.

###### pacs:

42.55.Px, 85.25.Cp, 74.45.+c, 42.60.MiAll numerous degrees of freedom of a bulk superconductor are frozen at low temperature. So its state is characterized by a single variable, the quantum phase of the superconducting order parameter.tinkham () In principle, this phase remains constant, that is coherent, for infinitely long time. Another example of a highly coherent system is a laser. There a macroscopic number of photons form a coherent state characterized by an optical phase and an amplitude.scully () In contrast to superconductors the coherence time is not infinite, but limited by the intrinsic noise in the device, caused by spontaneous photon emission. This leads to a phase drift and thus to coherence loss at typical time scales ofscullylamb () , with the number of photons in the mode and the cavity escape rate. Coupling these two very different systems may open up novel ways of controlling either of them, allowing for manipulation of one with the other.

Recently, a device has been proposedgodschalk () that combines the coherence of laser light with the coherence of the superconducting order parameter, dubbed the “half-Josephson laser” (HJL). This device was based on yet another recently proposed device, the “Josephson light emitting diode”recher:10 () (JoLED), which, it was shown, can be used for quantum manipulation purposes.hassler () The light in the JoLED is generated by a biased semiconductor p-n junction realized on a nanowire. To acquire a narrow emission spectrum of the junction, quantum dots are defined at the p-n interface. Using nanowires an efficient coupling can be achieved between semi- and superconductors.vanDam&Doh () By connecting both ends of the described p-n semiconductor nanowire to superconductors, superconductivity is induced at the quantum dot regions that alter the quantum dot states. Electrons and holes jump from the leads to the quantum dots and recombine there, emitting a photon. Even though the JoLED is biased with a voltage much larger than the superconducting gap, the superconducting correlations are retained because any direct, nonradiative charge transfer is blocked by the potential barriers in the p-n junction. The electron-hole recombination can happen in two ways. A “blue” photon at about the Josephson frequency, , is emitted when Cooper pairs from either side of the junction recombine. The phase of this blue photon is locked to the superconducting phase difference. A process resembling spontaneous parametric down-conversionspdc () allows for the emission of two “red” photons in stead of a single blue one, at about half the Josephson frequency, . The JoLED is the basis for the setup of the HJL. The latter exploits the emission of red photons by the JoLED and amplifies it by stimulated emission in an optical resonator. In contrast to the case of the JoLED, where red photons are emitted spontaneously, the red photons in the resonant mode of the HJL form a coherent state. The phase of this coherent state is locked to the superconducting phase difference.

Lasing was found for the red emission with the optical phase locked to the superconducting phase difference, where only two values of the optical phase are allowed, having a phase difference of . The importance of the phase locking in this device is that it can result in exponentially long laser coherence times. The idea to prolong coherence times by phase locking has been thoroughly explored in the field of lasers.ohtsu () Using, for instance, the technique of injection lockinglang () one could pump a laser with another one with longer coherence time, so that the pumped laser inherits this longer coherence time. However, in that case the coherence time is still limited by the fundamental noise processes in the pump laser, with some coherence time . The essential difference with the HJL is that the laser light inherits the coherence time of the superconductors, which is infinitely long. It will not become infinitely long, though, because the fundamental noise in the laser is still present and causes switchings to the lasing state with opposite phase, with exponentially small probability. Hence, the coherence times can be exponentially long.

Unfortunately, despite this promising prospect, the HJL does not exhibit these exponentially long coherence times because switchings between quantum dot states result in a shorter coherence time. The magnitude and phase of the laser light depends on the actual quantum state of the JoLED in the HJL. Not every state supports lasing however, and the ones that do give rise to various magnitudes and phases of the laser light. Switchings between the quantum dot states therefore cause decoherence of the HJL. The switching rates are suppressed because they occur in combination with the emission of off-resonant photons, which is allowed because of the finite linewidth of the optical resonator. For an optimal choice of parameters it is found that, despite this very different decoherence mechanism, the coherence time is, , similar to that of the common laser and much shorter than the envisioned exponentially long coherence times.

The design for the HJL is based on a single quantum emitter, putting rather high demands on the quality factor of the optical resonator. It was proposedgodschalk () to relax these demands by considering a large number, , emitters in a single resonator, increasing proportionally the upper bound on the damping rate at which lasing is achieved. As an added benefit, intensity fluctuations are expected to reduce as , thus realizing a less fluctuating light intensity.

In this article, we thoroughly investigate the alternative idea, where, instead of a single JoLED, we consider a large number of emitters forming a dipole moment to drive the optical resonator mode. The model that we introduce will be quite general as we will only consider macroscopic quantities and leave out the microscopic details. The crucial aspect of the design, the coupling of the quantum states to both superconducting leads, will of course remain, since it both drives the laser and induces the phase lock between optical phase and superconducting phase difference. The high demands on the quality factor of the resonator can be relaxed because it is driven by a big number of emitters.

We find that stationary lasing is possible with two equivalent values of the optical phase, as in the HJL, and with exponentially long coherence times owing to the phase lock. Rather surprisingly, the sole increase in number of emitters does not directly guarantee the increased number of photons in the mode. The average dipole moment of all emitters will turn out to be zero when the emitters are coupled to the resonator mode only. Incoherent emission to other modes is required to create a population imbalance,populationinversion () resulting in a nonzero average dipole moment. Owing to the large number of emitters, single switchings will not cause decoherence as in the HJL. Rather, they induce minor fluctuations of the dipole moment, transferring to detuning fluctuations in the resonator mode. These again transfer to fluctuations of magnitude and phase of the lasing state. The detuning fluctuations compete with the intrinsic fluctuations of the resonator mode. When both are sufficiently small, the superconducting phase difference will remain locked to a single value of the optical phase. The major decoherence mechanism of the laser is then provided by large fluctuations, represented by the tails of the noise distributions. These bring the optical field to the state with opposite phase. Large switchings occur rarely, on exponentially long time scales. The importance of such a feature is hard to underestimate. Indeed, examining these processes we have found coherence times (surpassing the fundamental limit for common lasers), scaling as , up to some critical value , where it saturates. While the linewidth of the laser is directly related to this coherence time, there are large tails in the line shape associated with the noise in the dipole moment as well as to the intrinsic quantum noise.

In the model we study, we neglect the noise source causing decoherence of the superconducting phase difference: voltage fluctuations. The superconducting phase difference is fixed only for ideal voltage bias. In case of nonideal bias, the voltage fluctuations will drive the phase evolution and cause the drift of the superconducting phase difference. The corresponding coherence time can be estimated using the voltage correlation function for a junction with impedance

(1) |

with the thermal energy for temperature . With this the average value of the phase is given by

(2) |

where the second Josephson relation and the voltage correlation function were used for the final expression. From this we find a coherence time of . The coherence time can thus be prolonged by reducing temperature or impedance of the nanostructure. Alternatively one can increase the number of Josephson junctions to reduce the noise. The latter approach is used in metrology to realize the Josephson voltage standard. A voltage drift of a few nV per hour can be realizedvoltagestandard () in that way, leading to coherence times on the order of s. Hence the voltage fluctuations may be disregarded when studying the decoherence in the device.

This article is organized as follows. In Sec. I we derive the governing equations for a broader class of devices, characterized by a large number of emitters. We start by briefly treating the model of the HJL, as described in Ref. godschalk, , and subsequently argue that the resulting equations can be generalized to represent a wide range of devices. In Sec. II it is shown that lasing with a phase lock occurs for these general equations in certain parameter regimes. There are distinguished regimes, separated by first and second order phase transitions. In Sec. III conditions are derived to determine the parameter regime where the noise can be considered in linear approximation. The frequency dependence of the noise, when in a stationary lasing state, is calculated. The noise in linear approximation cannot cause the decoherence of the laser owing to the phase lock. However, occasional large fluctuations will cause decoherence. These are treated in Sec. IV, for which the formalism of optimal paths is used. Finally, we conclude in Sec. V.

## I Model

The purpose of this section is to derive a rather general but simple model describing a novel class of lasing devices, those are driven by superconducting leads biased at voltage , have the optical resonance frequency at half the Josephson generation frequency, , and contain many quantum states capable of emitting light. The purpose of the general approach is to capture the essential properties of this class of devices, without involving specific microscopic details of the used model. In this section we start by briefly repeating the essentials of the model for the HJL at microscopic level and extend it to contain a large number of effective emitters, yet at the microscopic level. We formulate the generic features of the model which are not specific to the details of the microscopic description. Finally, we present the general but simple model where all these generic features are incorporated. We formulate the model in terms of a Fokker-Planck equation for three variables.

### i.1 The HJL with a single emitter

The HJL is described in Ref. godschalk, as an optical resonator driven by a single quantum emitter. The emitter is a double quantum dot in a p-n semiconductor nanowire, capable of containing up to two electrons and two heavy holes. The energy of the corresponding discrete quantum states is determined by orbital and Coulomb energies and is represented by Hamiltonian . Owing to the superconducting proximity effect, a gap is induced at the quantum dots for both electrons and holes, represented by Hamiltonian

(3) |

where the are the proximity induced pair potentials, with the phases retaining the values of the corresponding superconducting leads. Owing to gauge invariance, only the superconducting phase difference, will have a physical significance. Interaction of the double quantum dot with the resonator mode occurs according to

(4) |

where is the frequency detuning from resonance, with the resonance frequency of the optical resonator and the bias over the device. The operator , taken in rotating frame, represents the electron hole recombination required for the photon emission. It can be seen as the dipole moment operator that drives the resonator mode. Coefficient is the coupling strength of the dipole moment to the resonator mode. It is given by the dipole strength owing to quantum fluctuations of the electric field in the resonator mode. The interaction Hamiltonian is assumed to be perturbative as to allow for a large photon number before nonlinearities become manifest, which will induce saturation of the dipole moment.

The dynamics of the resonator mode can be treated in a semiclassical fashion, because a large number of photons will occupy the mode. In view of this, the , and will be replaced by their respective expectation values in the photon-dependent part of the full Hamiltonian. A scaled quantity is defined, . Within this approximation, the expectation value of the dipole moment, , is found by taking the derivative of the spectrum to , with and eigenstates associated with the full Hamiltonian, . Hence, the value of the dipole moment, , depends on the particular eigenstate of the double quantum dot. On the one hand, the dipole moment depends on the optical field, , but on the other hand, the evolution of is determined by the dipole moment. This leads to a set of self-consistent equations

(5) |

Stationary lasing is found when for finite values of . This state must be stable against small perturbations.

To get lasing, the rate at which photons are created in the resonator mode, , needs to exceed the cavity escape rate, . In common lasers, the imaginary part of the susceptibility, Im, results from pumping to excited states. Usually, to achieve lasing, the pumping must be sufficiently intense to invert the population. For the HJL, it is the superconductivity that creates the required complex susceptibility so that the population inversion is not necessary. In fact, the ability of the superconductors to absorb or supply Cooper pairs undermines the notions of groundstate and excited state in the quantum dots. For instance, the induced pair potentials enable the possibility of transitions from lower to higher energy states, while emitting a photon, with the required energy provided by the bias. This feature allows for cycles of emission or absorption of photons, in order to (de)populate the resonator mode.

Decoherence in the HJL with a single emitter is caused by spontaneous switchings between quantum dot eigenstates. Owing to the finite cavity linewidth, emission of off-resonant photons occurs inducing these switchings They occur at a slow rate, , and therefore introduce extra dynamics in the HJL described by a master equation. This is in addition to the already present dynamics of Eq. (5). It was found in Ref. godschalk, that even for an optimal choice of parameters, not every quantum dot eigenstate couples to the resonator mode strongly enough to support lasing. The ones that do support lasing will yield different stationary values . Since only some states support lasing, the optical field may even extinguish after a switch, to turn on again after another one. Hence, the switchings change , causing large intensity and phase fluctuations. Most notably, for every lasing state there is another one at . When going from some lasing state to a nonlasing state and returning to the same lasing state again, the optical field can go to both with equal probability. Hence the switchings are responsible for decoherence. The coherence time is determined by the switching rate, , and is similar to that of the common laser.scullylamb ()

### i.2 The HJL with many emitters

Even though the HJL, as described above, has a rather specific design, the framework provided by Eqs. (3)-(5) is quite general. It applies to a class of devices where a biased Josephson junction contains a structure of which the eigenstates couple to both superconducting leads and which can emit light by electron hole recombination. The light is emitted into an optical resonator. When the structure in the Josephson junction is large enough it can be divided into blocks. Every block represents an effective quantum emitter, which is a quantum system with a number of discrete eigenstates that can emit light. Every emitter must couple to both superconducting leads and must emit light into the resonator mode. In principle, the emitters interact with each other, but we may model them to be independent. As an example of such a structure, one might consider an ensemble of p-n semiconductor nanowires with a double quantum dot, the nanowire that was used in the previous section, or one might consider a two-dimensional electron gas between the superconducting electrodes.S-2DEG-S ()

Increasing the number of emitters in the Josephson junction to enhance the lasing is a method reminiscent of a line of research where superconducting lasers are made based on stacked Josephson junctions.stackedjj () There, the Josephson junctions are connected in series so that cascade amplification of the radiation occurs. The setup of the HJL is different though, because the emitters can be regarded as Josephson junctions connected in parallel rather than in series.

The dipole moment corresponding to the structure in the Josephson junction is determined by its energy, which is taken in a rotating frame of reference. When the structure contains emitters the total energy is given by the sum of the energies of the emitters. Assuming that emitter is in a state with energy , the total energy is simply . This yields a value of the total dipole moment given by , which should be inserted in Eq. (5) to look for solutions that support lasing, .

Every emitter is subjected to spontaneous switchings that change its eigenstate. Therefore, at time scales longer than the switching times it is not the total energy (dipole moment) that is relevant but the average total energy (dipole moment). The switching dynamics defines a master equation for each of the emitters, where the probability of emitter to be in state evolves as

(6) |

with the transition rate from state in emitter . The quantity in brackets defines a transition matrix, which necessarily has a null vector to ensure the conservation of probability. The null vector gives the stationary occupation probabilities, . With these, the average total energy is given by , where is the energy of state in emitter . Because of the coupling to the resonator mode, the energy depends on and the total average dipole moment is . Now this quantity determines the lasing solutions of Eq. (5).

### i.3 The average value of the dipole moment

Let us discuss a complication that might easily arise when considering this model or designing a device. In our approach we consider spontaneous switchings only to be a result of off-resonant photon emission. In the limit of a large number of emitters the associated transition rate matrix, Eq. (6), turns out to be approximately symmetric, with . The transition rates are calculated using Fermi’s golden rule. They are, up to some constants, given by the product of a matrix element squared, , for a given emitter , and a Lorentzian shaped density of states, , with the energy difference between the eigenstates and . We note that the photon number generically remains the same during a transition; is constant. Only after the transition may relaxation to the new stationary state occur. For a large number of emitters, however, a single switching is not enough to change the electrical field in the mode significantly. The stationary state remains almost the same after the transition and therefore also the eigenstates and energies remain approximately the same. Consequently, the transition rates also remain the same. Since the energies have a large mismatch with the resonance frequency (the emitted photons are far off-resonant) they do not depend significantly on detuning, , and we have , for fixed . The transition rate from state to is therefore equal to the transition rate of the reversed process. Hence, because the transition rates remain approximately the same after a transition, the transition rate matrix is approximately symmetric.

Innocent as it may seem, this feature of the transition matrix poses a problem for our device, because as a result the average value of the dipole moment, , must be approximately zero. The null vector for a symmetric transition matrix must be a vector with all entries equal. Therefore, in the stationary solution to the master equations, all the states in an effective emitter have the same probability to be occupied. The total average dipole moment is then given by . However, for the summation over we find

(7) |

where Tr is the trace over the space of states of emitter . The equality to zero follows because Tr is independent of . Consequently . We note that this is generic for any kind of system governed by our Hamiltonian and where the relevant states defining the dipoles, are equally populated. Naively one would expect that a large number of emitters would result in a large number of photons in the resonator mode. We see that this is not so. On the contrary, the total dipole moment could be zero and no photons are present.

We need to adjust the model used in Ref. godschalk, to get a nonzero average dipole moment. This is done by making the transition rate matrix asymmetric. Then the occupation probabilities for the states of an emitter will not be equal anymore, resulting in an imbalance in populations of states. There is no a priori reason to favor, for instance, the state with higher energy over the state with lower energy, for higher occupation probability. The notion of higher or lower energy of states has become obscured by the superconducting pair potentials. Hence, the created imbalance is not the same as a population inversion. Importantly, recombination of electron-hole pairs will still be characterized by the phase lock, because they are still coherent with the superconductors. Despite introducing incoherent transition rates, the special characteristics of the device are thus maintained.

An asymmetry in the transition rate matrix can come about in various ways. The simplest way is to take into account the photon emission into the environment, rather than to the resonant mode. If the rate of this incoherent process is comparable with the emission rate to the resonant mode, the average dipole moment is of the order of its maximum value per emitter. In a realistic system we expect relaxation transitions between the emitter states that do not involve photon emission. These relaxation channels will generally yield asymmetric transition rate matrices. One can also engineer transition rates to get a sufficiently asymmetric transition matrix. A straightforward approach is to change the rates that are already present by changing the electromagnetic environment. A switching to a state with higher (lower) energy is accompanied by emission of an off-resonant photon, with energy (), where . By changing the density of states in favor of one of these two, the spontaneous emission rate for the corresponding process increases with respect to the other one, thus creating the desired imbalance in transition rates. This might be done by inserting the device in a low quality factor resonator with resonance frequency about either of the frequencies. The low quality factor will ensure the rapid escape of the photon, before it is reabsorbed again. Alternatively, the electromagnetic environment can be changed by shining with a red or blue detuned laser on the device, favoring transitions to lower and higher energy eigenstates respectively, proportional to the laser intensity. This resembles the technique of resolved sideband cooling/heating.diedrich () Another possibility is to add new transition rates to the system by adding nonsuperconducting leads to the device and directly pump electrons and/or holes into the optically active region. This, however, leads to big challenges in the technical device design.

### i.4 Phase dependence of the energy

We have discussed the general features of possible microscopic realizations of a HJL with a large number of emitters. Let us construct a model that does not rely on microscopic details. The relevant quantity to deal with is the -dependent part of the total energy of the device, which determines the dipole moment. We show that the energy has a generic dependence on the phase difference between the superconducting leads.

Generally, we note that both common lasers and superconductors are systems characterized by a spontaneously broken -symmetry, a phase symmetry. As such they are characterized by an order parameter which has a magnitude and a phase. Any actual realization of the order parameter bears no apparent connection to the original unbroken symmetry. However, as a remnant of the symmetry, every value of the phase is equally probable and corresponding to states with equal energy. States with different phases are degenerate so that the energy of the system is independent of the phase.

Our model combines both of these systems, lasers and superconductors, which leads to a phase lock between the optical phase and the superconducting phase difference. Indeed, owing to the interaction between these two systems a phase dependence must enter the energy of the system so that not both of the spontaneously broken symmetries can be maintained. The emission of a single photon is caused by the recombination of an electron-hole pair, with the latter carrying the phase of their respective condensates. The phase of the terms in the relevant interaction Hamiltonian, Eq. (4), is therefore given by , with the optical phase and the superconducting phase difference. In fact, the full interaction depends on the cosine of this phase combination, leading to a similar dependence in the energy. As a consequence an energy minimum is found for a specific value of , which is a manifestation of the phase lock between optical phase and superconducting phase difference. A further inspection of the phase dependence of the energy shows that there must, in fact, be two equivalent values of the phases. Changing leaves the system invariant since the superconducting phase difference has not changed, but the phase in the interaction has changed to , so we conclude that . Hence, the original spontaneously broken symmetries of the Hamiltonian have been reduced to a single one and a spontaneously broken symmetry under rotations of of the optical phase. The remaining symmetry comprises invariance under a simultaneous change in optical phase and superconducting phase difference such that the quantity is unchanged. We stress that this reasoning applies to any system that combines two order parameters with spontaneously broken symmetries, with an interaction as in Eq. (4). Any such system manifests a phase lock.

### i.5 Energy

To find lasing in the general model, we need to have an expression for the energy of the device. This allows us to extract essential physics from the model.

Generally, the energy should depend on the value of light field in the resonator, , and the pair potentials, , in the superconductors: . We use parameter limits where the four variables are relatively small, enabling us to use a Taylor expansion for the energy. This limit corresponds to a weak drive of the resonator mode. To define “small” we go back to the model where the dipole medium consists of many emitters. We can associate an energy scale, , to each of these emitters. The induced pair potentials should be small compared to this scale, . In a realistic device this can be engineered by making suitable (tunnel) barriers between the superconductors and emitters. For the second limit, it is reasonable to take the interaction between dipole and resonator mode in the weak coupling limit. Hence the coupling constant of the interaction Hamiltonian, Eq. 4, is much smaller than . Therefore, this interaction perturbs the dipole energy by a value of the order of only when a large number of photons, , occupies the resonator mode. We want to make sure that the stationary value of photon number, , is always much smaller than , so that the use of an expansion of the energy in terms of is always justified.

When writing the Taylor expansion of the energy, , only a limited number of terms are allowed because of the specific phase dependence that is required. As was discussed in the previous section, the energy should depend on the phase combination, . Therefore, the lowest order phase-dependent term in the energy must be proportional to complex conjugate. This is a fourth-order term in the Taylor expansion. All other terms of equal or lower order are proportional to powers of or . Writing only the terms with immediate relevance, the energy is, up to fourth-order in Taylor expansion

(8) |

In this expression the phase of the induced pair potentials was absorbed into the optical phase, so it is redefined as . Furthermore, we absorbed the induced pair potentials in the prefactor , with a constant. The primes in and refer, respectively, to first and second derivatives at of the function , which is independent of the induced pair potentials. The lowest order term, , is a sum of and some other terms depending on but not on and . These terms will therefore not be relevant for the dipole moment. There are also terms proportional to , but these are small compared to , because they are of higher order in the Taylor expansion, and can safely be neglected. It is necessary to keep the term in the expansion, since this nonlinear term stabilizes a possible lasing instability. Without this term the light field would grow exponentially under unstable conditions. When calculating the dipole moment later, we will see that adds up to the detuning from resonance. Hence, the value of can implicitly be changed by changing the detuning. In fact, to get lasing the detuning should be such that effectively becomes small, of the order of or smaller, or goes to zero. Because the detuning , this also sets the scale . We stress that this is an effective change of . In the expansion of Eq. (8) the term is still considered to be larger than the fourth-order terms. With effectively , the term is of the same order as, or smaller than, the term so the latter should be maintained in the expression for the energy. With the above expression for the energy we have managed to separate the dependence of the energy on and from the microscopic details of the driving medium. The latter are all contained in the prefactors of the Taylor expansion.

Substituting Eq. (8) into Eq. (5), we find an explicit, general self-consistency equation. Noting that and , we have

(9) |

where and is the redefined detuning. This equation describes the deterministic evolution of the electric field, , in the resonator mode. We note that an equivalent complex conjugate equation exists.

### i.6 Fokker-Planck equation

The energy was considered to be a fixed quantity in the above section, corresponding to the average value of the dipole moment. To fully account for all dynamics in the HJL with many emitters, we also have to consider noise. This can be incorporated in a single Fokker-Planck equation that encompasses two noise sources. The first noise source is the traditional quantum noise due to photon emission from the mode.risken () The second noise source describes the fluctuations of parameters in Eq. (8).

To understand the latter let us note that Eq. (8) is valid not only in average, but also for each specific configuration of the emitter quantum states. The parameters , , and are configuration dependent. Since the emitters undergo spontaneous transitions the parameters fluctuate accordingly. Explicitly taking into account all the details of the fluctuations would greatly complicate the model. Instead, we employ the central limit theorem stating that any ensemble of equal distributions can be described as a Gaussian distribution when the number of distributions goes to infinity. Therefore, since we have a large number of emitters, as described above, we may regard the parameters as stochastical quantities with Gaussian distributions.

In principle, the noise in all parameters is expected to be of the same order of magnitude relative to their average values. However, since is of lower order in the Taylor expansion of Eq. (8), its fluctuations are more important than that of and . We therefore only consider a single noise source by taking and to be equal to their average values, while , with being the average value and being the deviation. The latter appears in Eq. (9) as an addition to the detuning. The noise can therefore be seen as fluctuations in detuning or in the resonant frequency. The corresponding variance, , is an important parameter of our general model. We shall take into account that emitters switch on a certain time scale that determines autocorrelations of these detuning fluctuations. The autocorrelation function is thus given byvankampen ()

(10) |

where is the inverse autocorrelation time of the emitter ensemble

It is possible to capture all relevant dynamics compactly in a single nonlinear Fokker-Planck equation. This is an equation for the probability distribution function of three real values, two of them are represented in the optical field and the third in the detuning deviation. The drift and diffusion terms in the equation respectively describe deterministic and noise driven, indeterministic processes.vankampen () The Fokker-Planck equation for our device encompasses two Gaussian noise sources, describing fluctuations in three quantities. First, there is the dynamics of which is subjected to a noise source, , describing the fundamental fluctuations in . This is the quantum noise associated with the discrete changes of photon number in the resonator mode. It satisfies, and . Second, as discussed, the noise in is equivalent to detuning fluctuations. This is represented by the quantity , with average and variance . The Fokker-Planck equation reads

(11) |

The single derivative terms on the right-hand side of Eq. (11) are the drift terms. Evaluating by using Eq. (11) will yield the self-consistency equation (9) for . We remark that in this context the brackets, , denote the average over the classical distribution , and not the expectation value of the operator . The diffusion is described by the second derivative terms on the right-hand side of Eq. (11).

Let us note that all the parameters we use have the dimension of frequency. It is convenient to make them dimensionless, measuring them in relative units of . This amounts to setting in Eq. (11). For clarity, we rescale the detuning deviation in Eq. (11) to appropriate units, by taking , which is the scale of the frequency domain where lasing is possible, as we show in Sec. II. Additionally we have . As a consequence, the term in Eq. (11) proportional to changes to . Other terms retain their form.

We conclude this section by stressing that this Fokker-Planck equation applies to a whole class of lasers, where a biased Josephson junction contains a structure of which the eigenstates couple to both superconducting leads and which can emit light by electron-hole recombination into the resonator mode.

## Ii Lasing

In the previous section we have derived a model for a class of lasers that is driven by superconductivity and that exhibits a phase lock. Here we show that lasing with a phase lock is indeed possible and we give expressions for the lasing threshold and the number of photons in the resonator mode. We find three different regimes, with a different number of stationary lasing states.

A full description of the dynamics of lasing would require solving Fokker-Planck equation (11), but this is difficult, if not impossible, even for a relatively simple system such as ours. Instead, we first search for stationary lasing solutions to the deterministic equation of motion, Eq. (9), thus disregarding noise. Once these solutions are found we study the effects of noise when near a stable stationary lasing state. This is the topic of the subsequent two sections.

To find the desired stationary solutions to Eq. (9), with , we rewrite them as real equations. Defining leads to a set of equations,

(12) |

where is the photon number in the resonator mode. These equations can be readily solved to find an expression for and for the phase of , which was defined as ,

(13) |

where we conveniently rescaled the detuning, . Hence, we find stationary lasing solutions. The fixed value of implies that and are locked. Note that the equations are invariant under the transformation and , implying the occurrence of two solutions which differ in phase by , as predicted in the previous section. Indeed, Eqs. (13) are invariant under the transformation . Hence and represent four stationary solutions. Besides these, there is another stationary solution, being , delivering a total of five possible stationary solutions. We remark that only depends on and not on and . Furthermore, a large photon number can be achieved when .

Studying the stationary solutions of Eq. (13) we find three different regimes that differ in the number of physical solutions. Since represents the number of photons in the resonator mode, it must be either a real and positive quantity or it must be zero. As a result we note that lasing solutions, with , can only be achieved when . Then, depending on , three different regimes exist (Fig. 1). We choose but similar regimes can be found when . For and , both so that the only physical solution is . We call the parameter regime with a single solution regime I. Regime II is given by the condition that , where and . Including there are now three physical solutions. The third regime is when and , where both , corresponding to a total of five physical solutions.

Even though there can be multiple physical solutions, not all of them will support stationary lasing, because some will be unstable against small perturbations. To find the condition for stability we expand Eq. (9), in terms of and , up to first order about a stationary value. This yields

(14) |

where and are at the stationary values, and are the respective deviations from the stationary values and Eq. (12) was used to rewrite some terms. When an eigenvalue of the matrix is negative (positive), a perturbation in the direction of the corresponding eigenvector decays (grows) exponentially to (from) the stationary value of and . The case where both eigenvalues are negative corresponds to a stable stationary lasing state. After a fluctuation in any direction, the system will relax to the stationary state again. With one positive and one negative eigenvalue there is a stable and an unstable direction against perturbations, thus defining a saddle point. The case where both eigenvalues are positive cannot occur, because the trace of the matrix is smaller than zero. The condition of the stability therefore depends on the sign of the determinant of the matrix, positive being the stable stationary state. Evaluating the determinant and again using Eq. (12) we derive that the stability for a physical solution is determined by sign, whereas the stability of the solution at is determined by sign.

To illustrate how the system relaxes to a stable stationary lasing state after a small perturbation, we give the eigenvalues of the matrix in Eq. (14). These are, for the case of , and, for the case of , . Hence, when is sufficiently large, the system oscillates about the stationary point with a frequency and relaxes to it in rotations.

For the stability of the stationary solutions in the various regimes is as depicted in Fig. 1. Regime II is the most relevant one for lasing as the solution is unstable, whereas the solution is stable. In regime III the and solutions are all stable, whereas the solutions are unstable. When , the qualitative picture is the same. Quantitatively, the location of the regimes changes. In the phase diagram of Fig. 1, the regimes I and III are interchanged so that the plots for are mirror images of the plots for . Additionally, the () solution is now the stable (unstable) one.

The boundaries between the regimes can be seen as lasing thresholds because they separate the lasing regimes from the nonlasing regime. The boundaries between regimes I and II and between regimes II and III are lines in parameter space where the determinant of the matrix of Eq. (14) is zero at . They are given by . For negative the value () corresponds to the boundary between regimes I and II (II and III). The boundary between regimes I and III is determined by the condition . At this boundary, the solution is still stable, but two pairs of stable and unstable solutions arise at values . When gradually crossing the boundaries from I to II or II to III, the system undergoes a second order phase transition [arrows (i) and (ii) in Fig. 1]. Thereby the state becomes unstable. In contrast to this, when crossing the boundary from regime I to III [arrow (iii) in Fig. 1] the system will undergo a first-order phase transition, because all physical stationary solutions will arise directly at nonzero values, two stable and two unstable. We note that cannot be increased indefinitely as the pairing potentials must remain small enough to satisfy the adopted approximations.

To conclude this section, we find that regime II is the only regime where lasing occurs without an interruption caused by a large fluctuation. This is because only there the solution is unstable. In contrast to this, in regime III a large fluctuation can bring the device from a lasing state to the stable nonlasing state at thereby interrupting the emission of light. Regime II occurs in a small interval at large detuning, where . The number of photons in the resonator mode, being of the order of , is large when . In the remainder of this work, we assume we have chosen the parameters and such that the system is in regime II. We take so that the solution is the lasing solution.

## Iii Noise and its frequency dependence

Let us consider fluctuations in the device. It is important to note that because of the phase lock the laser cannot lose coherence by a random drift of the optical phase as in common lasers. On the contrary, after a small fluctuation a damping force returns the system to the stable stationary lasing state. The average value of the optical phase is therefore retained. As in the HJL,godschalk () the only possibility to lose coherence is when the device switches to another stable stationary lasing state. This requires a sufficiently large fluctuation in photon number or optical phase. We envision these large fluctuations to occur only rarely, on exponentially long time scales, implying that the noise in must be linear, i.e. dominated by small fluctuations about the stationary lasing state. In this section we investigate the HJL in the linear noise approximation. There are two noise sources causing intrinsic quantum fluctuations in the resonator mode and fluctuations in detuning. Therefore, the steady state of the HJL can be described by a multivariate Gaussian probability distribution, for which we calculate the variances. Several parameter limits will be discussed concisely for the photon number variance. We also find an expression for the corresponding noise spectrum. The large fluctuations that cause decoherence of the laser are the topic of the next section.

We first describe an intuitive picture of the shape of the full noise spectrum, containing the spectrum of both large and small fluctuations. A sketch of this is given in Fig. 2. The large fluctuations occur at exponentially long time scales, , satisfying . Here, as we see later in this section, the small fluctuations occur at much shorter time scales, determined by the largest of the response time of the resonator mode, , or the autocorrelation time of the detuning noise, , yielding . Because of this separation of time scales, both kinds of fluctuation contribute independently to the noise spectrum. This is given by the Fourier transform of the correlator of ,

(15) |

At time scales the correlator is expected to decay because of the switching to another lasing state. The switching changes the optical phase but retains the stationary photon number . The correlator decays as

The noise spectrum of this correlator is a Lorentzian with exponentially small width, and with peak height . At time scales small fluctuations have a dominating contribution to the correlator. Writing for a small fluctuation, , we have . The requirement that the fluctuations are small implies that they are Gaussian distributed as shown below. One expects the noise spectrum to consist of three peaks. One peak is associated with the fluctuations in detuning, and has the width . Two side peaks, shifted with respect to zero detuning, are associated with the intrinsic quantum fluctuations and have width . These peaks are shifted because of the oscillatory behavior of small deviations from equilibrium position, which was described in the previous section in the context of stability analysis. In fact, the shifts should be . The noise spectrum approximately corresponds to a sum of three Lorentzians. Here the peak heights must be . Adding up the two contributions, the total noise spectrum is a high and very narrow peak on a background formed by three low and broad peaks, as depicted in Fig. 2. In this section we concentrate on the small fluctuations that are responsible for the broad background in the noise spectrum of the figure.

A full description of the small fluctuations can be found by linearizing the Fokker-Planck equation of Eq. (11). This amounts to expanding to and only retaining up to linear terms in , , and . The resulting Fokker-Planck equation can be solved exactly,vankampen () but we do not do this, because we are only interested in the steady state probability distribution. The latter is fully characterized by the variances of the variables of the system.

Instead of directly using the linearized Fokker-Planck equation of Eq. (11) to calculate the variances, we use the equivalent system of Langevin equations for a set of appropriate variables. The advantage is that in the course of this procedure, we directly obtain the noise spectrum associated with the linear noise. The linearized form of Eq. (11) is fully equivalentvankampen () to the following system of Langevin equations,

(16) | ||||

Here it is noted that the first equation is complex and represents two equations. There are two different Langevin forces satisfying

(17) | ||||||

with other correlators being zero. The source represents the quantum noise in the resonator mode, while the source represents the noise in the detuning. The current set of equations is not the most convenient one because the variables and are not orthogonal to each other. It will be most natural to work in a system with a variable , representing radial displacements or photon number changes, and a variable , representing azimuthal displacements or phase changes. The relation between these two bases comes about in the following way. One can write the deviation in vector form in a Cartesian coordinate system as . Alternatively, in a cylindrical coordinate system , where the follows because is taken about , which is separated from the origin with a distance . The variations represent the variations in optical phase, but the variations correspond to the fluctuations of the square root of . Hence, . Using the transformation relations, and , one can derive expressions for and to find that

Using this we can rewrite the above set of Langevin equations. Reminding the reader that we use for the lasing solution and [Eq. (13)], we arrive at

(18) | ||||

Here we have three Langevin forces satisfying

(19) | ||||||

The other correlators are zero.

With the change of variables, we now need to find the steady-state probability distribution to the Fokker-Planck equation for the variables , , and . This is given byvankampen ()

(20) | ||||

where is simply the variance of the distribution of the detuning fluctuations. Now using the Langevin equations we can find the variance matrix according to the following procedure. First, we take the Fourier transform of Eqs. (III). For any variable we have . The resulting system of equations is solved for , and . To get the variances, we calculate the inverse Fourier transform at of the expectation value of products of the Fourier transformed variables. For instance, we have

(21) |

According to this procedure, we find the variances and write them in a form that will be convenient later:

(22) | ||||

We remark that some variances diverge when approaching the lasing threshold, which is when or . Also, the combination is equal to the determinant of the stability matrix of Eq. (14). A larger value of the determinant implies higher stability, as the variances decrease.

With the expressions for the variances, we are able to quantify in which parameter regimes the linear Fokker-Planck equation and the steady-state distribution of Eq. (20) are valid. Their validity is ensured when and , so that all expressions in Eq. (III) must be much smaller than one. Note that under these conditions, there is no a priori bound on , which is already purely linear. It will be particularly useful to study the expression for in a bit more detail, because it is always of the order of, or larger than the other ones. Additionally, this is an interesting quantity, because large fluctuations in lead to switchings and are thus important for the long coherence times.

The variance of the fluctuations in photon number is studied in a few limits. For reasons of simplicity and clarity, it is most convenient to be sufficiently far from the lasing threshold, so we take and . Keeping this in mind, the first limit to be studied is the one where . Then we have

(23) |

We note a few properties. To start with, the term proportional to is associated with the noise in the detuning, whereas the term is associated with the intrinsic quantum noise of the resonator mode. The latter exactly corresponds to the fundamental noise expected in any common laser.scully () In line with this division, there are two distinct regimes when considering the dependence on . For small enough values of the intrinsic quantum noise dominates so that the variance decreases inversely proportional to . Increasing beyond a critical value the variance saturates as the term becomes unimportant and the noise in detuning starts to dominate. We further note that the magnitude of the fluctuations is independent of and independent of in the limit . The criterion to have small fluctuations now leads to the requirements that and .

The two previously mentioned limits of in the detuning noise dominated regime can be intuitively understood. If , it is smaller than the decay rate of the resonator mode so that the fluctuations in are slow enough for the optical field to adiabatically follow them. The number fluctuations are then purely determined by the equilibrium dipole distribution. When , the optical field cannot follow the dipole fluctuations adiabatically. However, because , the optical field will respond fast to the dipole fluctuation before it fades away. In fact, the optical field will start to rotate about the stationary point with a period and the fluctuation only fades away after rotations about the stationary lasing state. The variance of the number fluctuations is now determined by the instantaneous, out-of-equilibrium value of the fluctuation of , as determined by in Eq. (17), leading to an extra factor of compared to the previous case with .

We continue by studying the variance of photon number in the second limit, which is and leads to

(24) |

In this limit, again there is in the dependence on a division in two regimes: one where the variance is dominated by the noise in detuning and one where the intrinsic quantum noise of the resonator mode is most important. In contrast to previous case there is now a dependence on both and . Here, the criterion to have small fluctuations leads to the requirements that and .

There is also an intuitive way to understand the detuning noise-dominated regime in the limit . Here, the dynamics of the detuning fluctuations is on much shorter time scales, , than the dynamics of the resonator mode, , so that the latter has a reduced time window to respond to the detuning fluctuations. Compared to the earlier case with , the reduction factor is .

As part of the calculation done to find the variances, we can also find an expression for the noise spectrum, which is defined as . Defining we have

(25) |

The term proportional to is directly related to the noise in the detuning, whereas the other term is directly related to the quantum noise of the resonator mode. Which of the two terms is the most important depends, of course, on the parameter regime, as discussed above. The noise spectrum can be written in three parts, as we show in the Appendix and discussed in the intuitive picture above, with each part corresponding to a peak. There is a central Lorentzian peak at with width and there are two non-Lorentzian sideband peaks at . In the limit and the sideband peaks become Lorentzian with width . In this limit, the central peak typically dominates over the sideband peaks when . For , the central peak almost vanishes compared to the sideband peaks.

We have now described the linear fluctuations in the HJL and derived conditions on parameters for which fluctuations are small, . The noise in various parameter regimes was studied. For small , but still much larger than , the intrinsic quantum noise of the resonator mode was found to have a dominating contribution to the photon number fluctuations, while for large , the noise in the detuning is the most important. There, the magnitude of number fluctuations depends on and the ratio . The noise spectrum corresponding to the small fluctuations was found and has a central peak with varying width, and two sideband peaks with width .

## Iv Large Fluctuations

In the previous section we studied small fluctuations in the HJL under the assumption that they dominate the fluctuation spectrum. These small fluctuations do not cause decoherence of the laser. Decoherence is caused by large fluctuations that switch the stationary lasing state to one with a phase difference of ; see fig 1. In this section we study the large fluctuations thoroughly to find the time scales at which the switchings would occur. This gives the decoherence time of the HJL.

The probability to have a large fluctuation can be calculated using the formalism of optimal paths.reviewop () This formalism is, for instance, used in classical diffusion driven systems where the transition rate from one stationary state to another, across a high barrier, is calculated. The path that crosses the barrier and has the largest probability is called the optimal path. Its shape and probability can be found using the principle of least action. Its probability determines the transition rate.

In this section, we introduce the formalism of optimal paths and apply it to our system, after which some generic properties of the action are derived. After this we show how the coherence time is related to the optimal paths. We investigate the dependence of the action of the optimal paths on the device parameters. This is done both by using the results of previous section and by simulating the optimal paths.

### iv.1 Trajectories and relation to Kramers’ escape problem

To start with, we outline the dynamics of the device in the absence of noise, describe the effect of noise and explain the relation to Kramers’ escape problem.

In Sec. II we have derived the requirements on stable stationary lasing in the HJL. Let us suppose that the optical field in the HJL is not in a point of stable lasing and that . If we disregard the noises we know that the field will evolve to one of the points. This evolution is described by Eq. (9). We visualize the dynamics by making a stream plot, as is done in Fig. 3. To describe the dynamics it is convenient to associate the field components, , with coordinates of a “particle” subject to a coordinate-dependent “force field” that causes the motion of the particle with the velocity proportional to the force. Starting from a given initial condition, the particle will flow towards a stable stationary point. Its trajectory corresponds to a streamline in the plot.

Figure 3 corresponds to regime II as defined in Sec. II. In the figure the three stationary points are indicated by , , and . The point , at is the unstable saddle point, whereas the points and are stable points or attractors. For a general initial condition in the vicinity of the saddle point the particle is repelled from this vicinity. There is, however, a particular line called stable direction such that the particle is attracted to the saddle point if the initial condition is chosen to be on this line. The streamlines coming to the saddle point at the stable direction form the separatrix of Eq. (9), which is the boundary separating the domains of attraction of the attractors and . We see that the particle circles around an attractor before reaching it. This signifies that the force field cannot be represented as a gradient of a scalar potential. The force field must also have a transversal component that causes the rotation and can be represented as the curl of a vector potential.

Full equations of motion shall also include a variation of detuning . With this we have three coupled equations of motion that define the force field in the space of three variables,

(26) |

The last equation gives the relaxation of the detuning to its equilibrium position. The specifics of this equation is the presence of a saddle line. The line is a saddle line because for each value of the point is a saddle point with respect to the motion in directions.

In Sec. III noise was introduced. With the Langevin noise added to the dynamics, Fig. 3 no longer gives a full description of the particle evolution. The particle in the ‘force field’ will experience random kicks, causing it to change from one trajectory to another, even when initially at rest at an attractor. Effectively, the kicks cause the particle to diffuse away from its starting point, even against the streamlines of the force field. Nonetheless, the force field counters the diffusive flow trying to bring the particle back to the attractor. These two kinds of motion have a fundamentally different nature. The evolution along the flow lines is deterministic, in contrast to the diffusive motion which is nondeterministic. Importantly, the diffusive motion of the particle lifts up the restriction of the particle to stay in the domain of attraction of an attractor. A series of random kicks may repel the particle from an attractor and eventually pull it over the separatrix into the domain of attraction of the other attractor. After that the particle relaxes naturally along a streamline to the other attractor. This trajectory corresponds to the large fluctuation we are after: A switch between the two stable lasing states has occurred.

Let us suppose, in contrast to what is the case for Eq. (9), that a particle moves in a potential well, of which the gradient defines the corresponding force field, and that the fluctuation-dissipation theorem is applicable, like in Kramers’ escape problem. In this case the Fokker-Planck equation can be solved to calculate the escape probability from the potential well.kramers (); vankampen () This escape probability involves the Arrhenius factor, , where is the height of the barrier that is crossed by the particle, counted from the potential minimum, and is the effective temperature of the particle. However, because of the transversal component of the force field associated with Eq. (9), we cannot use this general approach to solve Fokker-Planck equation (11). There, the escape probability does not depend on a potential barrier height, but on the shape and probability of the path that crosses the separatrix. We note that the potential barrier and the separatrix both form the boundary of the domain of attraction where the particle is trapped.

### iv.2 Optimal paths and the principle of least action

With the conceptual framework in mind, let us continue by showing how the probability to have a transition from one stationary state to the other, is related to the concept of the optimal path. We only consider paths that start in stationary point and end in stationary point . The time elapsing during the transition will be infinite, going from to , because the particle can only leave and approach asymptotically.

It is possible to express the probability distribution for the fluctuations in terms of an action. We start with a general set of Langevin equations, written in vector form,

(27) |

The Langevin forces represent Gaussian white noise, satisfying and , thus supposing an infinitesimally short autocorrelation time. The noise can be expressed as a single probability distribution that is a product of the probability distributions of the separate Langevin forces at each moment in time,

(28) |

with the normalization constant. One may interpret this expression in another perspective, as the probability distribution for histories of fluctuations, because the integral in the exponent represents the history of fluctuations. Each history corresponds to a unique progress in time of . We stress that the values of at each time remain strictly independent from the ones at other times . So how does this concept of histories aid us? Even though the Langevin forces might have infinitesimally short autocorrelation times, the vector field will need some finite time to relax to the stationary state. During this time, the history of fluctuations will be relevant. The resulting cumulative effect of fluctuations is best captured by using Eq. (27) to rewrite the exponent of , yielding

(29) |

Hence, we indeed see that the noise probability distribution can be described using the action, defined as .

The probability distribution, , is strongly reminiscent of a path integral. Indeed, is a probability density which, upon integration over all paths, gives the probability to go from to . When fluctuations are small, as is described in Sec. III, the path integral is dominated by the “classical” path and by the quadratic fluctuations about it.kleinert () The classical path is found by the principle of least action. It is what we have called the optimal path. The optimal path corresponds to the lowest value of the action for the transition, , which relates to the transition probability

with the prefactor corresponding to the quadratic fluctuations. Hence, the probability of a transition from to depends on the optimal path.

Without specifying any details about we note two things. First, by inspection we see that , with the equality satisfied for the deterministic path, for which . Second, the integrand of is equivalent to a Lagrangian and because it only implicitly depends on time, a conserved quantity similar to the Hamiltonian exists

(30) |

The equality of the Hamiltonian to zero is not universal, but is appropriate for the path we consider. This path connects the stationary points, where and therefore . Being a conserved quantity, must always be zero at the path connecting the stationary points.

The action as it is given in Eq. (29) is not in the most convenient form. To minimize this action we need to find both the shape of the optimal path and the speed with which the path is traveled. By exploiting the conserved Hamiltonian, we can relate the speed to the force field resulting in an action where the optimal path is only determined by its shape. To show this, we choose the diffusion constants, without losing generality, as for all . This is allowed, because we can always rescale and accordingly. Applying this to Eq. (30) we then find that . With this, the terms and can be rewritten and added to each other to become . Because both terms in the Lagrangian now only contain a single time derivative, we can parametrize time at will. We parametrize it as a monotonously growing function , such that and . Consequently , so that the action becomes

(31) |

Minimizing this action will result in the optimal path irrespective of the exact form of parametrization and therefore irrespective of the speed with which the optimal path is traveled.

From the expressions for the action, Eq. (31), an intuitive understanding of the dynamics of the system can be obtained. Suppose