Theory of the coherence of topological lasers

We present a theoretical study of the temporal and spatial coherence properties of a topological laser device built by including saturable gain on the edge sites of a Harper--Hofstadter lattice for photons. In the comoving frame of its chiral motion, the lasing edge mode behaves very similarly to a one-dimensional system, and the correlations of long-wavelength fluctuations display the peculiar Kardar-Parisi-Zhang (KPZ) scaling. At longer times, the finite size of the device starts to matter, and the functional form of the coherence decay changes from a KPZ stretched exponential to a pure exponential. Still, the nonlinear dynamics of KPZ fluctuations reflects into an enhanced linewidth as compared to a single-mode Schawlow-Townes phase diffusion. While all these features are general to extended laser arrays, the role of topology in protecting the coherence from static disorder is finally highlighted. This opens exciting possibilities and intriguing open questions both for fundamental studies of non-equilibrium statistical mechanics and for concrete applications to laser devices.

The laser is one of the most fundamental tools in modern science [1,2]. Its defining feature is the emission of radiation characterized by an unprecedented long coherence length and time [3]. This makes laser sources essential ingredients in a wide range of applications, which justifies a continuous theoretical and technological research of new devices. Also froam a fundamental science perspectives, the laser mechanism represents also one of the archetypical models of nonlinear physics, non-equilibrium statistical mechanics, and quantum optics [4][5][6][7].
Recent years have witnessed an explosion of the field of topological photonics. Following the ground-breaking of a number of work [8,9], experiments have proved the chiral propagation of modes along the edges of passive photonic latticesdisplaying topologically non-trivial band structure in different material platforms [10], and then their self-sustained lasing, first for 1D chains of polariton resonators [11] and then in 2D photonic crystal [12] and microcavity arrays [13]. While the robustness of the lasing operation in the presence of a strong static defect on the edge or moderate on-site disorder has already been highlighted in the above experiments and theoretically characterized in [14,15], a dynamical study of the temporal and spatial coherence properties of the topological laser is still missing, and it will be the topic of this work.
In the absence of external noise sources, the long time coherence of a single-mode laser is ruled by the diffusion of the phase of the macroscopic electromagnetic field, so that the corresponding correlation function features an exponential decay. As it was first theoretically understood by Schawlow and Townes [3], the diffusion of the phase and the corresponding emission linewidth of an ideal device are induced by spontaneous emission events. While this intrinsic noise mechanism represents the ultimate quantum limit for the coherence of the laser operation, the diffusion process itself can be treated by a semiclassical approach based on stochastic differential equations [5,16]. Given the importance of this Schawlow-Townes result in standard laser theory, it seems natural to undertake a similar analysis for a topological laser device. This is the main subject of the present work. For the sake of simplicity, we focus on the simplest scenario of a so-called class-A laser [17], where we can neglect the dynamics of the carriers and restrict the model to a coherent laser field amplified by an instantaneous gain medium with linear refractive index, and subject to a white noise.
Since lasing occurs into a chiral mode localized on the one-dimensional edge of the lattice, we expect that many properties may be similar to 1D arrays of coupled laser resonators, for which it is known [18][19][20][21][22][23] that, for small intensity fluctuations, the long-wavelength phase dynamics follows a Kardar-Parisi-Zhang (KPZ) equation [24]. Here, we numerically prove that the same holds for the chiral mode of the Harper-Hofstadter topological laser model studied in [15], once the phase fluctuations are studied in the reference frame moving at the group velocity of the chiral mode. In particular, the spatio-temporal correlation functions are shown to follow the universal KPZ scaling function with no dramatic renormalization of the nonlinear coupling.
While the KPZ scaling holds either in spatially infinite systems or at intermediate times, the long-time behaviour of finite systems is determined by the spatial extension of the lattice. In the textbook theory of lasers, the Schawlow-Townes line is computed by considering the phase diffusion of a single mode [16,25] or, from a different but equivalent perspective, by looking at the drift induced by noise on the Goldstone collective mode of the linearized Bogoliubov theory [26,27]; also, a nonlinear refractive index and non-homogeneity of the lasing mode can be straightforwardly included [28]. However, while in both approaches the coherence time turns out to be proportional to the length of the system, one has to keep in mind that they are neglecting the KPZ nonlinearity. This can in fact give a further contribution to the linewidth by generation of effective noise from finite momentum excitations into the zero-momentum global phase mode. While a similar physics was considered at the Beliaev-Landau perturbative level in theoretical studies of the phase diffusion of atomic Bose-Einstein condensates [29], to our knowledge no such calculation has ever been performed for the linewidth of a long array of coupled laser resonators. Our numerical calculations highlight a crossover for growing system size from a single-mode Schawlow-Townes scaling of the coherence time, proportional to the total photon number, to a reinforced diffusion due to the nonlinear dynamics of spatial fluctuations.
The structure of the article is the following. In Section I we focus on finite 1D laser chains and, after inspecting the crossover from KPZ to diffusion at very long times, we illustrate the significant reduction of the coherence time from the single-mode Schawlow-Townes prediction due to nonlinear spatial fluctuations of the phase. In Section II we first review the basic concepts of the topological laser device based on the chiral edge mode of a 2D Harper-Hofstadter lattice. We then investigate its emission spectrum and the one-dimensional character of its phase fluctuations. In particular, we probe the universal KPZ dynamics of long wavelength fluctuations and we compute the linewidth of the laser emission at very long times. In Section III we demonstrate the robustness of topological laser operation and of its coherence in the presence of static on-site disorder. Conclusions and perspectives are finally drawn in Section IV.

I. 1D ARRAY OF LASER RESONATORS
In this Section we review the coherence of extended lasers, also called driven-dissipative quasicondensates [6], in the semi-classical approach. While a lot of work has focused on the connection with the KPZ equation and its infrared limit [18][19][20][21][22][23], not much attention has been paid to the equal-space time correlation function for finite systems, one of the most relevant quantities for experiments. Since in what follows we will be massively concerned with this quantity, we devote some effort to its study.
We start by presenting a model for a single laser resonator and by recalling the well-known Schawlow-Townes result for the linewidth. We then move to consider finite arrays of coupled resonators.

A. Single resonator linewidth
The basic features of laser operation inside of a resonator, namely the breaking of the U (1) symmetry of the field and stabilization via saturation, are caught by the semi-classical Ginzburg-Landau equation where ψ is the macroscopic electromagnetic field in the resonator, the intensity of the field in the resonator is n = |ψ| 2 , γ represents the losses of the resonator, P the gain and n S a saturation coefficient. The external forcing ξ is white noise ξ * (t)ξ(t ) = δ(t − t ) with diffusion coefficient D. In the absence of the external noise Eq. (1) has the U (1) symmetry for the phase of the field; the steady-state field ψ 0 is zero below the lasing threshold P < P th = γ, while in the symmetry broken lasing phase it satisfies n 0 = n S P P th − 1 . In particular, for P = 2P th it holds the very transparent expression n 0 = n S .
To deal analytically with Eq. (1) it is convenient to resort to the modulus-phase formalism: writing the field as ψ = √ ne iφ one has where now one has two real uncorrelated noises ξ l (t)ξ l (t ) = δ l,l δ(t − t ), l, l = 1, 2. For small enough perturbations, the intensity relaxes to n 0 with a Lyapunov rate Γ = γ(P −γ) P ; on the other hand, as a direct consequence of the U (1) symmetry of the model, the dynamics of the phase is a diffusion with no restoring force.
Assuming noise is small enough to cause minor perturbations to the intensity and assuming phase differences are gaussianly distributed (cumulant approximation [18]), the autocorrelation of the field reads so that the usual Brownian motion scaling entails that the decay of coherence is described by the exponential In Fourier space, this corresponds to a lorentzian power spectral density. The coherence decay rate is the celebrated Schawlow-Townes linewidth [3,25]. The crucial feature of this formula is the intensity appearing at the denominator, connected with the following picture: the more photons in the resonator, the less the phase of the field is perturbed when a photon with a random phase is spontaneously emitted into the field.

B. Extended lasers and KPZ equation
For a single resonator, the linewidth is inversely proportional to the number of photons. For an array of coupled resonators, the naive intuition is that the number of photons within some correlation length of the system will play a role. In particular, for small system sizes all resonators oscillate in phase and, for a given number of photons per resonator n 0 , the linewidth scales as the inverse of the length of the array. Instead, for large enough systems and finite observation times, it is known that the physics is described by the Kardar-Parisi-Zhang (KPZ) equation [21].
We start by generalizing Eq. (1) to a one-dimensional array of N x coupled laser resonators x = 1, ..., N x with coupling J and ξ * x (t)ξ x (t ) = δ xx δ(t−t ). Periodic boundary conditions are also assumed. For J large and calling a the distance between neighbouring resonators, one can interpret Eq. (7) as a discrete version of a continuous field of mass m given by J = 1 2ma 2 . The corresponding continuous equation is the complex Ginzburg-Landau equation [26]; assuming fast relaxation of the intensity fluctuations, it turns out that the dynamics for the phase is described by the Kuramoto-Sivashinsky equation (KSE) where φ here is the unwound phase living on the real axis, and not the compact one restricted to [0, 2π].
The characteristic scales of the system as a function of the microscopic parameters have been reported in [18]: Measuring space, time and (unwound) phase in terms of l * , t * , φ * the adimensional KSE reads Since for the 2D topological laser to be discussed below other scales can be relevant, in most of the paper we will use the physical (without tilde) space, time and phase. A rescaling will be proposed in order to study KPZ features in Sec. II D.
The renormalization group analysis shows that at long distances and times the KSE flows to the KPZ universality class [30]. The KPZ equation [24] reads and it was originally proposed to describe the growth of interfaces. Its scaling behaviour at low energies (and assuming an infinite system and stationary regime) is characterized by the correlation function and by two exponents χ, z which determine the asymptotic behaviour of the spatial and temporal correlations respectively, according to ∆φ 2 x,0 ∼x 2χ and ∆φ 2 0,t ∼ t 2χ/z . In 1D we have χ = 1/2 for the roughness exponent and z = 3/2 for the dynamical one; even more precisely it holds where A = D/ν is the variance of ∂xφ. The universal function g KP Z is known exactly [31] and we just recall that g KP Z (u) − 2|u| → 0 for u → ∞ and g KP Z (u) → 1.150... for u → 0. As a consequence, the equal-time correlation function has the random walk form ∆φ 2 x,0 = A 2 |x|, for which the KPZ nonlinearity λ doesn't play any role. In other words, only looking at the spatial correlations the dynamics is not distinguishable from the linear Edwards-Wilkinson one, with χ = 1/2 in both cases; instead a difference will be recognizable in the temporal correlations, with z = 2 in the linear theory and z = 3/2 in KPZ. In fact, at small distances and times, the nonlinear term in the KSE can be neglected and the linear, or Bogoliubov, analysis applies. It has also been shown [20] that the crossover from Edwards-Wilkinson to KPZ physics is slower in the presence of a nonlinear refractive index.

C. Linewidth of extended lasers
In a photonic experiment one has easily access to the equal-space time correlator where we dropped the dependence on x in assuming a uniform system. For a finite system and very long observation times, the coherence decay has to be at least as fast as a pure exponential g (1) (t) ∼ e −|t|/τc , with τ c the coherence time. This is a consequence of the U (1) symmetry, hence the diffusion of the phase. As we are going to review, the Schawlow-Townes coherence time holds in the Bogoliubov limit, where the linearization on top of the noiseless state is a good approximation and in particular the phase variation along the system can be considered small [27]. Within linear response modes of different momenta are decoupled and for a discrete lattice (16) whereD k are the effective drifts coefficients, connected to the shape of the Bogoliubov modes, and tend to D for k → 0. At small momenta, the Bogoliubov frequency of the diffusive branch reads ω k −iγ k = −iJ 2 Γ −1 k 4 [26]. Then, for long times only the k = 0 mode contributes as ∼ |t|, since the other modes decay exponentially, giving corrections of order γ −1 k |t|. The factor 1/N x represents the fact that white spatial noise is equivalent to drawing noise realizations with a given k and the probability to pick k = 0 is one over the number of modes. The same conclusion holds also for finite continuous systems by replacing N x with the length of the chain (i.e. the total number of photons n 0 N x matters). On the other hand, for infinite (continuous or discrete) systems the sum is replaced by an integral to yield the Bogoliubov prediction g (1) (t) ∼ e −B|t| 3/4 with B a constant; the slow decay comes from the k = 0 noise mode occurring with probability zero. In the presence of a nonlinear refractive index one has to use γ k ∝ k 2 , leading to the Bogoliubov prediction g (1) (t) ∼ e −B|t| 1/2 for the infinite continuous wire.
Of course, we know that the Bogoliubov approximation is not adequate for an infinite or large enough system, where KPZ features set in. For an infinite system the stretched exponential behaviour g (1) (t) ∼ e −B|t| 2β can be proven [21], with 2β = 2χ/z = 2/3 and B nonuniversal. If the system is large but finite, the Bogoliubov approximation breaks down, but U (1) symmetry still ensures the very long time coherence to decay at least as fast as a pure exponential, g (1) (t) ∼ e −|t|/τc . The KPZ physics typical of the infinite chain is still visible at an intermediate timescale, up to a saturation time scaling as N z x . To illustrate these arguments, in Fig. 1 we report the time correlation function computed by solving numerically Eq. (7) while monitoring a single site of a laser chain, for 3 different sizes N x = 64, 256, 1024, panels (ac) respectively. More precisely, the thick black line is − log g (1) (t), which measures the phase drift, and the red and cyan lines are linear fits in the loglog scale of the plot. Keeping the same observation window, for small sizes the decay is mainly diffusive, while for larger sizes the Schawlow-Townes behaviour is pushed at very long times and on a finite time scale the KPZ stretched exponent g (1) (t) ∼ e −B|t| 2/3 is seen. In order to find in a very clean way the exponents of both regimes, one would need very large sizes and arbitrary observation times, which is numerically too demanding. Here we use a small hopping J = 0.5γ to enhance KPZ physics and P = 2γ, n s = 1000, D = γ; simulations are initiated with the uniform field ψ 0 = √ n 0 .
While linear theory Eq. (16) decorrelates different kmodes and predicts the Bogoliubov-Schawlow-Townes coherence time Eq. (15) proportional to the total number of photons n 0 N x (or n 0 L x for the continuum), the nonlinearity of the KSE plays a role for large enough systems. This is seen in the numerical results for the coherence time τ c of chains of increasing sizes and is reported in Fig. 2 divided by the single site Schawlow-Townes coherence time τ ST,1 = 2n0 D . The blue and cyan points, corresponding to coupling strengths J = 5γ and J = 0.2γ respectively, follow the Bogoliubov result until a certain critical size (that from numerics seems to scale with J as l * ∼ J 4/7 ), after which they deviate with a very peculiar behaviour. A deviation is indeed expected and it can be understood looking at the KSE or KPZ equation: the total phase drift can be decomposed in two uncorrelated contributions, one generated by the k = 0 component of noise that yields the Schawlow-Townes drift, and one coming from the phase field components at finite k, which scatter via the nonlinearity λ into an additional k = 0 drift. Still, the observed non-monotonicity of τ c (N x ) suggests that a non-perturbative analysis is needed for a quantitative explanation of the phenomenon, that will be the object of futue work. A perturbative study of the phase drift of equilibrium condensates involving Beliaev decay was presented in [29].
Notice that the required numerical effort grows up rapidly with N x , since one has to study very long times, larger than the KPZ saturation time ∼ N 3/2 x . For this reason, a statistical analysis of the errors has been restricted to N x = 128, J = 5γ, where we performed 8 runs of duration T = 10 8 γ −1 . For each run we extracted τ c ; the average of these τ c 's is very close to the coherence time of the averaged g (1) (t). More specifically, the standard deviation of the single run τ c is comparable with the size of the marker in Fig. 2. For all other N x we perform a single run of the same duration T ; consequently, the variability of the single shot is smaller for smaller N x 's, and larger for N x = 256, the largest size investigated here.

II. HARPER-HOFSTADTER TOPOLOGICAL LASER
As shown in recent theoretical [14,15] and experimental [12,13] works on photonic platforms, it is possible to make the edge mode of a topological insulator to lase by introducing a gain, maintaining at the same time topo- logical properties such as the chirality of propagation and robustness in circumventing defects. We devote the rest of the paper to the semi-classical study of phase fluctuations of the field in the presence of external white noise, accounting for quantum and classical noise. Studying the chiral laser operation in the reference frame of the edge mode, we find strong analogies with the physics of 1D laser arrays reviewed in the previous section. In particular, a discussion of the effective couplings characterizing the KPZ regime is provided.

A. Model
As in [15], the underlying topological insulator is provided by a 2D Harper-Hofstadter model. By labeling with x = 1, .., N x and y = 1, .., N y the lattice sites, the equations of motion for the field read i∂ t ψ x,y = −J ψ x,y+1 + ψ x,y−1 + e −2πiαy ψ x−1,y + +e 2πiαy ψ x+1,y + i 2 P δ y,1 1 + n x,y /n S − γ ψ x,y + 2D x,y ξ x,y , where α is the flux per plaquette of the synthetic magnetic field per plaquette in units of the magnetic flux quantum. Only α = 1/4 will be studied here. In this work we will consider a cylindrical lattice with periodic boundary conditions along the x direction, and we introduce the gain medium on the sites of the y = 1 edge. The diffusion coefficient of the noise is taken uniform D = γ; notice that the minimal amount of noise for a lossy resonator is set by the unitarity of quantum mechanics and corresponds to a semi-classical drift of D = γ/2 [5]. A space dependency of the diffusion coefficient, as it will be used later, does not change qualitatively the results.
Let' s now recall the properties of the Harper-Hofstadter Hamiltonian, momentarily neglecting gain, losses and noise. Within the Landau gauge used in Eq. (17), implementing periodic boundary conditions along the x direction makes the system translationally invariant, so that it is convenient to label states by the wavevector k x . For the lattice on a cylinder, α = 1/4 corresponds to having 4 bands consisting in states delocalized all over the bulk of the system, plus two chiral modes per side localized on the edge. On the y = 1 side per example, one mode has positive group velocity and negative eigenenergy, and viceversa for the other mode. The fact that there are no two counter-propagating states on the same edge at the same frequency is at the origin of the topological protection. This eigenstructure is depicted in Fig. 3.a, cyan-dotted lines.
Starting the simulation with zero field, noise breaks the U (1) symmetry and a positive or negative frequency mode is randomly selected for lasing. This physics has been investigated in [15] in the absence of noise but random initial conditions. The lasing modes are selected with a probability distribution depending on the discrete frequencies corresponding to the eigenvalues of the finite non-interacting system. In terms of the free Harper-Hofstadter spectrum, this distribution is symmetric with respect to zero frequency and has support in the bulk gaps of the topological insulator. Its peak corresponds to the eigenstates that are most localized on the edge, since the gain is more effective for these. As pointed out in [32], the overlap of the Harper-Hofstadter eigenmode with the amplifying region of the lattice determines a k x -dependent imaginary part of the dispersion for the bare lattice. This remains true also for the Bogoliubov edge modes of the linearized theory on top of the lasing state [33].
The main features of the topological laser operation are illustrated in Fig. 3. In the lower panel we plot the field modulus |ψ x,y | for a N x = 32, N y = 12 cylindrical lattice,showing localization of the mode on the edge. The upper panel reports the power spectral density S(k x , ω) of the field on the edge ψ x,1 (t), from which one can see the two modes, with opposite chirality, living on the y = 1 side. The lasing mode, denoted by the cyan dots, dominates the spectrum, with the N y = 12 bare Harper-Hofstadter bands in dotted lines. The brightness of the gap modes comes from localization on the edge, where the medium is amplifying; a full account based on Bogoliubov analysis of all these features will be presented soon [33].

B. Chiral fluctuations
In this paragraph we study the phase fluctuations of the 1D non-equilibrium condensate defined by the field living on the amplifying boundary of the Harper-Hofstadter lattice, ψ(x, t) ≡ ψ x,1 (t). In particular, for visualizing the physics it is convenient to make a change of reference frame to be comoving with the lasing chiral edge mode. We stress that in this paragraph we want to provide an illustrative picture in terms of the phase of the field, but all results in the rest of the paper are based just on the correlation functions of the field.
As seen above, lasing occurs in the edge modes living in the spectral gap of the bulk. The dispersion ω em (k x ) of this branch has a curvature J ef f (k las x ) < J and is characterized by a group velocity v g = dωem dkx (k las x ). We can determine k las x by counting the winding number of the phase; the average frequency ω las with which the phase of the field revolves can be extracted by fitting the field evolution of single sites. These spatial and temporal windings are somehow uninteresting, since they are not due to fluctuations of the field but to the deterministic dynamics of the device. To remove these and concentrate on stochastic fluctuations, we define the slowly varying field ψ sl (x, t) ≡ e −ik las x x+iω las t ψ(x, t). Looking at the phase of a typical realization of ψ sl (x, t), as shown in Fig. 4.a, a clear pattern can be seen, where a phase fluctuation moves with velocity close to v g and gets distorted. To visualize the intrinsic dynamics of how these phase bumps develop, we plot in Fig. 4.c the phase of the field in the comoving frame For these parameters (the same as the previous paragraph) with J = 5γ and N x = 128, one can notice that the phase fluctuations develop very slowly and are a quite small. To increase their relative importance one can turn on the noise (or equivalently turn down the mean intensity n 0 ), consider longer lattices or decrease J, as it will be discussed in Section II D.
While the transformation to ψ CM (x, t) allows for a direct visualization of the phase dynamics, it is possible to study the fluctuations circulating along the edge by computing the space-time correlation function of the original field, where the average is over noise realizations and invariance under time and x translations is assumed, and no estimate of k las x , ω las is required. As apparentlooking at the smooth stripes in Fig. 4, analysing g (1) (x, t) is the cleanest way to obtain the velocity at which fluctuations travel, which is perfectly compatible with v g .

C. Frequency comb
The spectrum |ψ(k x , ω)| 2 (where Fourier transform is implied) shown in Fig. 3.b can be understood by inspecting the Bogoliubov modes on top of the noiseless lasing state. This will be fully discussed in a work in preparation [33], but we can anticipate that, as a consequence of the U (1) symmetry breaking, the diffusive Goldstone branch, analogous to the one in [26], plays a major role. In Fig. 5.a we show the spectral density

S(ω) = 1
Nx x |ψ(x, ω + ω las )| 2 for the parameters and lasing point shown in Fig. 3.b. Apart for the main peak, a comb structure with symmetric side-peaks is visible. It turns out that the frequency spacing is determined by the discretization of the HH spectrum, i.e. v g 2π Nx . The visibility of the comb is not merely due to the existence of the states at those energy which are populated by noise, but it is enhanced by correlations: in Fig. 5.b we report the normalized momentum-space intensity-intensity correlation function where the averages here are over noise and time and momentum-space densities are evaluated at the same time over all the edge, n kx (t) = |ψ(k x , t)| 2 . For the first two side-peaks one can see nearly perfect correlation, indicating that the peaks are populated in pairs by parametric scattering processes; therefore, the sidepeaks are more intense for smaller values of the curvature of the Harper-Hofstadter dispersion J ef f at the lasing point. For uncorrelated excitations k x + k x = 2k las x one has R (2) (k x , k x ) = g (2) (k x )g (2) (k x ) −1/2 , with as usual g (2) (k x ) = 1 for k x = k las x and g (2) (k x ) = 2 otherwise. Similar generation of phase-matched photons occurs p.e. at the inflection point of the lower exciton-polariton branch [34]. Also the topologically trivial 1D chain of Section I displays similar combs when the lasing point is not k = 0 but close to the inflection of the cosinusoidal dispersion.

D. KPZ space-time correlations
The KPZ universal dynamics occurs, in a lattice of N x sites, on timescales shorter than the saturation time ∼ N z x , after which the Schawlow-Townes like behaviour described above sets in, but larger than the timescales where linear Bogoliubov and non-universal physics dominate. This requires the system to be large enough, precisely it should be at least N x √ 2πl * [18]. Then we consider a system of length N x = 1024 (N y = 12), which requires a small hopping J = 0.5γ in order to clearly observe the KPZ exponents, while keeping intensity fluctuations within 15%. Notice that the hopping J, hence the topological gap, is comparable with the bare singlesite linewidth γ, so that it is a priori not obvious that chiral edge modes may exist. Still, above threshold the laser linewidth is not directly set by γ and lasing in the edge mode is observed in the simulations; we have also verified the robustness of chiral lasing in the presence of a strong defect on the edge.
To measure the KPZ features we have performed 30 simulations of duration γT = 5 · 10 5 , starting from the wavevector for which the Harper-Hofstadter eigenstate is most localized on y = 1, that is k x = −2π 155 1024 . For each run the correlation function ψ * (x, t)ψ(0, 0) is computed and then averaged over the 30 trials to yield g (1) (x, t). The typical dynamics occurring in a time γT = 2 · 10 4 is depicted in Fig. 6.a, where the phase of the field ψ CM (x, t) along the edge and in the comoving frame is shown; the fractal structure typical of interface growth can be recognized.
Defining the correlation function in the comoving frame g (1) CM (x, t) = | ψ * CM (x, t)ψ CM (0, 0) | and the rescaled correlator KPZ universality requires In particular, one hasC(0) = A. In Fig.6.c we plot C(x,t) for x = ±20, ..., ±130 and demonstrate that in excellent approximation they depend only ont/x z [35]. The orange and cyan series correspond respectively to the topological laser and to the 1D chain lasing in the k = 0 mode. In both cases we used J = 0.5γ, N x = 1024, n S = 1000, P x,y = 2γδ y,1 , D = γ 2 (1 + δ y,1 ). The crucial point about Fig.6.c is that phase, space and time have been rescaled by φ * , l * , t * obtained by using the effective masses and gain parameters: so for the topological laser J ef f is the curvature of the Harper-Hofstadter band at the lasing point and P ef f is chosen so that a 1D chain would reproduce the observed intensity on the edge. For instance, J ef f 3.19γ hence l * 1.92 here. A detailed report of possible mappings of the 2D topological laser to a 1D one will be presented soon [33].
As anticipated, the Renormalization Group analysis [30] predicts for the 1D array lasing in the k las = 0 mode that the rescaled KSE Eq. (11) flows to the low energy effective KPZ theory Eq. (10) with λ = 2, since, thanks to the Galilean invariance holding for KSE and KPZ equations, the nonlinear coupling is not renormalized. Indeed, in this case the rescaled correlation functions collapse to a unique curve, which is decently fitted in Fig.6.c (blue dashed line) by using λ = 2, A = 0.96, as expected from the Galilean invariance argument.
Also for the topological laser λ = 2, A = 1.14 yields a good fit (red dashed line); as upper bound for the fit of the nonlinearity we mention λ = 2.14, A = 1.10. Notice that for the topological laser there is no guarantee that a Galilean invariant KSE holds microscopically; on the contrary, an analysis on the lines of [18,19] suggests that, while rescaling with J ef f still yields a microscopic λ =  Fig. 1 for the topological laser. Thick black lines correspond to the logarithm of the temporal correlation for a given site on the edge of the lattice. Thin dashes represent − log g (1) CM (t), that is the behaviour of coherence in the frame comoving with the edge state. For increasing system sizes and a given temporal window, a crossover between diffusive (red fits) and KPZ (cyan fits) decay of the coherence is observed. The amplitude of the oscillations is inversely proportional to the spatial coherence of the device; for the long time coherence only the global phase matters and the oscillations fade away.
2, other terms of the kind ∇ 2 φ(∇φ) 2 , (∇ 2 φ) 2 , ∇ 3 φ∇φ should be added in Eq. (10). These additional terms come from effective imaginary derivatives due to the k x dependent localization of the Harper-Hofstadter eigenstates on the edge, and we would expect them to become important at very low momenta, contributing to a higher coherence of the device; in particular, they break Galilean invariance and potentially renormalize significantly λ, since as the KSE term they contain four derivatives. This renormalization, though, is not seen in the numerical results.
Finally, we give a brief remark on the experimental protocol to assess KPZ physics. The analysis of the correlations functions g (1) CM (x, t) of Fig.6.b is carried over in the frame comoving with the chiral mode, while in the lab one can easily measure correlation functions between different points at different times g (1) (x, t). In the end it holds g (1) CM (x, t) = g (1) (x − v g t, t), so that g (1) CM is obtained by a straightforward post-processing of g (1) (x, t); graphically, this amounts to fit v g from data similar to Fig. 4.b and tilt the correlation function so to have at any t the maximum of g (1) (x, t) in x = 0.

E. Temporal coherence and linewidth
To complete our analysis, we concentrate on the temporal coherence properties of the device at very long times, or equivalently we study with high resolution the spectral linewidth of the main peak of Fig. 5.a. Mimicking a possible photonic experiment, we monitor the field of only one resonator on the amplifying side of the lattice and repeat the analysis of Section I B. We remark that, even though the phase fluctuations travel chirally and consequently at short times the correlation function of a given resonator displays some oscillations like in Fig. 4.b, for the long time coherence these fade away and are irrelevant and only the leading exponential decay matters. Similarly, monitoring a resonator not on the edge one sees a fast initial decay, but at long times the phase coherence is the same as the edge, since it is the global phase of all the lattice to be U (1) symmetric and to diffuse. As seen for the 1D chain, the scaling of the linewidth with N for small enough N x is caught by the Schawlow-Townes line computed in the Bogoliubov approximation, while at higher N x nonlinear KPZ-like nonlinearities will determine a deviation from a straight line.
To start with, we illustrate once again the crossover between Schawlow-Townes and KPZ coherence decay, by reporting in Fig. 7 the phase drifts − log g (1) (t) and − log g (1) CM (t) ≡ − log g (1) CM (0, t) in a given temporal window and increasing system sizes N x = 64, 256, 1024. In particular, given the importance of temporal correlations in distinguishing between the linear Edwards-Wilkinson and KPZ regimes, we remark that the hallmark of the dynamic exponent z = 3/2 is already visible in the simple photonic experiment consisting in the measurement of the single site g (1) (t). The sharp local minima of − log g (1) (t) are local maxima of the coherence and occurs with a period equal to N x /v g , the time for a fluctuation to make a round trip of the edge; indeed, we can see that these minima sample the equal space coherence function g (1) CM (t). The agreement of Fig. 7.c with the KPZ result is very good. Notice also that the oscillations of − log g (1) (t) are present at small times and are large if the spatial coherence of the device is poor; for the long time coherence, instead, only the global phase of the field over all the lattice matters and no oscillations occur.
Focusing now on the very long time exponential decay, the diffusion rate of the phase in small topological lattices can be computed in the Bogoliubov approach. In full generality for a nonlinear system in whatever dimensionality on a lattice of N sites labelled as x, let us call L las the 2N × 2N Bogoliubov matrix of the linearized dynamics on top of the lasing steady state. Let V = {V xσ,p } be the invertible matrix which diagonalizes L las , where σ =↑, ↓ is a spin notation for the particle and hole components of the Bogoliubov problem and p labels the 2N eigenmodes. The Goldstone mode V xσ,G , that we assume to be unique with all other excitations having a finite life-time, is the eigenstate with zero eigenvalue. The effective noise acting on the lasing mode will be determined by the projection of the bare noise on the mode itself at each point; for generality we consider a position dependent bare noise D x . Then, in the Bogoliubov approximation, the phase drift associated with the Goldstone mode is independent of the location on the lattice and given by where the summation represents the projection of the noise on the Goldstone mode and is actually independent of x, since Φ G is the normalization of the phase variation of the Golstone mode which is constant. If V V † = 1 a very clear expression holds for the Schawlow-Townes line: with n tot = x n x . For a topologically trivial uniform system lasing at the bottom of the band, the k = 0 sector of L las is diagonalized by a 2×2 unitary matrix V V † = 1, since modulus and phase are decoupled at k = 0. In the topological case the k x = 0 sector has dimension 2N y × 2N y and V is not unitary V V † = 1. Still, for the parameters considered it turns out that V V † 1, so that the reader may keep in mind the approximate expression In standard laser theory the non-orthogonality of V is also known under the terminology of Petermann factor and excess noise [28,36]. The numerically observed coherence time for the topological laser is reported in in Fig. 2, red dots; in this case the theory line refers to Eq. (23), with the full matrix V , not the unitary approximation. Qualitatively, the topological laser behaves similarly to the 1D laser array. As apparent for the first two red points , the Bogoliubov phase diffusion coefficient matches the numerics for small N , while a marked deviation occurs already for N = 64.
Notice that, while for small systems the Schawlow-Townes trend of the linewidth does not depend on J ef f nor on the small k x Bogoliubov modes (only k x = 0 matters), the deviation from the Bogoliubov theory observed for large systems instead does strongly depend on J ef f , in a way related to the KPZ nonlinearities illustrated above.

III. LASING WITH ON-SITE DISORDER
While in the previous paragraphs a theory of the fluctuations of the 2D topological laser was developed by exploiting the many similarities with a non-topological 1D array of lasers, in this Section we show that the two models have a dramatically different behaviour when some static disorder is present in the resonators. In particular, while for the topologically trivial laser a small amount of disorder strongly affects coherence in a realization dependent manner, the temporal coherence of the topological laser is preserved in a robust way.

A. 1D lasing with disorder
Since a certain degree of fabrication imperfections and inhomogeneities are always expected in realistic photonic or polaritonic lattices, we start by considering the effect of on-site disorder on the lasing properties of a topologically trivial array of resonators. We are interested in understanding if topology qualitatively changes the response to .... Establishing the relationship of our model with the theories developed by the broad Random Laser [37] community goes beyond our present goals.
Along all the Section, a disordered potential is added to Eqs. (7) and (17), where G(0, 1) is a Gaussian random variable with mean 0 and variance 1. Here we restrict to N x = 128 and J = 5γ. The lasing dynamics in the presence of disorder is in general very complex, but, since our ultimate goal is a qualitative comparison with the topological laser, we focus here on the coherence time of the system for various values of disorder W , and in particular on whether there is a threshold after which coherence collapses. For linear waves or a quantum problem, the sensitivity of the eigenstates at a given energy to a static perturbation is proportional to the spectral density of states. Then, in order to have a more fair comparison, we consider both the 1D chain lasing in k las = 0 and in k las = π 5 , which has a group velocity (hence a density of states) comparable with the central part of the Harper-Hofstadter edge mode.
As shown in Fig. 8, a small disorder has certainly some impact on the coherence time of the device, here normalized to the value in the absence of disorder. A detailed description of all the possible behaviours goes beyond the scope of this work, we just notice that in most of the simulations coherence is quite soon lost, while in few cases we even have a significant enhancement of τ c . In particular, we fit τ c for each site (here we use an exponential fit even if the shape of g (1) (t) is in general very complex) and plot the average over the lattice. For strong disorder W ∼ 0.8γ different portions of the sample lase with different frequencies, like hinted at in Fig. 9. Here the phase of the field is plotted in time and space (in the time rotating frame correspondent to the bottom of the band with no disorder) for W = 0.8γ; the field was initialized in k las = 0 and let evolve for some time before making the movie.
Notice that each shape of the marks in Fig. 8 indicate that the same realization of G(0, 1) and the same initial conditions have been used, and we vary only W . It is then particularly curious the non-monotonicity of τ c with W .

B. Topological robustness
The same protocol was repeated for the topological laser with N x = 128, N y = 12, J = 5γ. Results are reported in red squares and triangles in Fig. 8. Also, the crosses refer to a N x = 52, N y = 14 lattice with open boundary conditions and gain along all the edge [15], which contains the same total number (128) of sites.
In contrast to the 1D chain, the behaviour of the topological laser when disorder is introduced is quite regular and similar between different realizations, with temporal and spatial coherence being preserved until a certain threshold around W ∼ J, that is of the order of the Harper-Hofstadter gap. For small W , disorder has a minor impact at all. For intermediate values there is a surprising and systematic enhancement of the temporal coherence; from the movies of the phase of the field (SM) it emerges that fluctuations, instead of moving around the edge with the group velocity of the chiral mode, stay localized. The coherence time is anyway less than the ultimate Bogoliubov-Schawlow-Townes bound (23) which for W = 0 is τ ST /τ c 2.7; a possible interpretation, to be further investigated, is that disorder hampers. the dynamics leading to the detachment of points in Fig. 2 from the linear approximation line. At higher disorder, spatial and temporal coherence are eventually lost.
Moreover, while in the previous paragraph we initialized the laser in the k las = 0 or in k las = π 5 modes, here initial conditions are the vacuum of the field, demonstrating that disorder do not affect the capability of the topological laser to reach by itself the operational regime. It is also important to notice that results are approximatively independent of the realization of disorder.
In conclusion, the correlation time of the 1D chain is basically unpredictable in the presence of modest lattice inhomogeneities, while the behaviour of the topological laser is quite regular and robust.
FIG. 8. Temporal coherence: topological robustness. Coherence time, normalized to its value for a clean sample, for a few realizations of on-site disorder in the 1D chain laser (blue marks for lasing point k las = 0, green for k las = π 5 ) and in the 2D topological laser (red). Same marker shapes correspond to the same disorder potential and various height W .

IV. CONCLUSIONS
In this work we have investigated the spatio-temporal coherence properties of arrays of coupled laser resonators, focusing on analogies and differences between a nontopological 1D chain and chiral edge-state lasing in a 2D topological Harper-Hofstadter lattice.
In the 1D case we have highlighted the crossover for growing observation times from a Kardar-Parisi-Zhang scaling to a phase-diffusion regime. Also, the phase diffusion rate at long times displays for a growing system size a crossover from a single-mode Schawlow-Townes linewidth to a faster decoherence determined by the nonlinear dynamics of spatial fluctuations. Provided one works in the reference frame moving at the group velocity of the chiral mode, these results are found to directly apply to the chiral laser emission in the edge states of a topological laser device. In practical experiments, all these features can be observed by measuring the different-time correlation functions between different sites of the lattice in the lab frame.
While for clean samples the spatio-temporal correlations behave very similarly for the topological and trivial devices, topological protection entails a much larger resilience to inhomogeneities. For the 1D chain, static disorder is able to spatially localize the lasing mode and/or break it into several disconnected and incoherent pieces; on the other hand, the chiral motion of the edge state of a topological laser device is able to maintain the spatial and temporal coherence across the whole sample up to stronger values of the disorder strength of the order of the topological gap.
These results open exciting perspectives both for technological applications and for studies of fundamental physics. On the application side, we extend the robustness of topological lasing in the presence of static disorder, already established as far as the slope efficiency is concerned [14], to the coherence properties.
From the theoretical point of view, ongoing research includes the classification of the different Kardar-Parisi-Zhang universality subclasses for our model [23], the extension of our study to class-B semiconductor lasers by explicitly including the carrier dynamics [17] and Kerr nonlinearities [33], and the generalization of the topological laser concept to photonic lattices with different dimensionalities and different band topologies [10].
From the experimental side, the effectively periodic boundary conditions for the chiral motion of the edge mode due to its circulation around the lattice are extremely promising to suppress undesired spatial inhomogeneities and boundary effects in experimental studies of the critical properties of different non-equilibrium statistical models [38].