# Effects of correlated and uncorrelated quenched disorder on nearest-neighbor
coupled lasers

Amit Pando  Department of Physics of Complex Systems, Weizmann Institute of
Science, Rehovot 7610001, Israel  Sagie Gadasi  Department of Physics of
Complex Systems, Weizmann Institute of Science, Rehovot 7610001, Israel  Eran
Bernstein  Department of Physics of Complex Systems, Weizmann Institute of
Science, Rehovot 7610001, Israel  Nikita Stroev  Department of Physics of
Complex Systems, Weizmann Institute of Science, Rehovot 7610001, Israel  Asher
Friesem  Department of Physics of Complex Systems, Weizmann Institute of
Science, Rehovot 7610001, Israel  Nir Davidson  Department of Physics of
Complex Systems, Weizmann Institute of Science, Rehovot 7610001, Israel

###### Abstract

Quenched disorder is commonly investigated in the context of many body systems
such as a varying magnetic field in interacting spin models, or frequency
variance of interacting oscillators. It is often difficult to study the effect
of disorder on these systems experimentally since it requires a method to
change its properties in a controlled fashion. In this work, we study the
effect of quenched disorder in the form of frequency detuning on a coupled
lasers array using a novel degenerate cavity with tunable disorder and
coupling strength. By controlling the properties of the disorder such as its
magnitude and spatial correlations, we measure the gradual decrease of phase
locking due to the effects of disorder and demonstrate that the effects of
disorder depend on the ratio between its correlation length and the size of
the phase locked cluster.

Introduction -  Many different physical systems, both quantum and classical,
are well described by many-body interacting oscillators: Transverse field spin
models describe spins effectively rotating around the local magnetic field,
which can synchronize to reach finite magnetization even in the presence of a
spatially varying magnetic field [1, 2, 3, 4, 5]. The synchronization of
classical phase oscillators has been theoretically studied for decades through
the Kuramoto model[6], and is manifested in many different systems such as
arrays of Josephson junctions[7, 8], phase locking of coupled laser arrays [9,
10], and even human networks [11]. In all of these, disorder in the system is
one of the main obstacles to synchronization, acting against the interaction
between the individual members of the ensemble.

While many theoretical studies of these systems investigate how they are
affected by disorder, it is often hard to do so in experiments, as it requires
disturbing the experimental system in a controlled and accurate fashion. In
this work, we use a highly controlled system to introduce quenched disorder to
an array of nearest-neighbor coupled lasers in order to study its effects on
phase locking. The disorder is introduced in the form of frequency detuning,
where the resonant frequency of each individual laser in the array is shifted.
By precisely controlling the magnitude of the disorder, we show how it
gradually decreases the ability of the lasers to phase lock. By varying
spatial properties of the disorder, we are able to demonstrate how its effects
depend on the interplay between the scales of the problem, namely the
correlation length of the disorder and the size of the phase locked cluster of
lasers.

In addition to experimentally measuring these effects, our work also
highlights the potential of our experimental system for studying disorder:
Coupled laser arrays have been shown to display many interesting physical
phenomena, such as topological states[12, 13, 14], non-hermitian dynamics[15,
16], and geometric frustration[17], and our system allows to study them in the
presence of precisely controlled disorder. Furthermore, these results have
practical implications for managing the effects of quenched disorder on
physical systems, as physical realizations of such disorder often have an
inherent correlation time or length relating to the process from which they
originate.

Experimental Arrangement -  To study the effects of detuning disorder on
coupled oscillators, we use a modified digital degenerate cavity laser
(DDCL)[18, 19, 20] to form 400 coupled lasers in a 20 by 20 square lattice.
The experimental system, shown schematically in Fig. 1, includes a 4f
telescope, an ND:YVO4 gain medium, a spatial light modulator (SLM), and a
tunable coupling mechanism. The gain medium is pumped by a quasi-CW $808nm$
diode laser with a pulse duration of $500\mu s$ and a 4Hz repetition rate. The
SLM is used as a digital phase and amplitude mask, allowing us to create a
laser array with arbitrary geometry and to control the frequency detuning
between lasers: By changing the phase on the SLM pixels, we vary the effective
cavity length and hence the resonant frequency of each laser with a precision
of $\Delta\Omega=\frac{\pi}{256\tau_{c}}$, where $\tau_{c}=\frac{2l}{c}\approx
13.3ns$ is the cavity roundtrip time. The diameter of each lattice site of the
amplitude mask is $d_{laser}=200\mu m$ and the distance between neighboring
lattice sites is $d_{lat}=300\mu m$.

Figure 1: Schematic drawing of the experimental arrangement. The intra-cavity
SLM is used to form the laser array via a square lattice amplitude mask and to
introduce detuning to each laser. The tunable coupling is highlighted by the
dashed line: By rotating the $\lambda/2$ waveplate, we can change the
magnitude of the coupling between the lasers continuously. The insets show
representative FF (left) and NF (right) intensity distributions. The four
sharp peaks indicate phase locking in an out-of-phase configuration ($\pi$
phase difference between neighbors) due to the negative sign of the coupling.

We introduce tunable negative nearest neighbor coupling between the lasers via
modified Talbot coupling[21, 22]: We insert a polarizing beam splitter (PBS)
into the cavity such that it deflects a portion of the propagating light
towards an identical amplitude mask located at the image plane of the
telescope, and an additional mirror which is displaced by an additional
$\frac{1}{4}Z_{talbot}=\frac{d_{lat}^{2}}{2\lambda}$ from a 4f telescope
configuration. The rotation angle $\theta$ of a $\lambda/2$-waveplate inserted
into the cavity determines the amount of light deflected to the coupling
branch of the cavity, effectively changing the coupling magnitude without
changing its other properties. We define the coupling as a normalized overlap
integral[23, 24]:

$K_{mn}=\frac{\innerproduct{E_{m}^{i-1}}{E_{n}^{i}}}{\innerproduct{E_{m}^{i-1}}{E_{m}^{i}}},$
(1)

where $m$ and $n$ are the laser indices, and $i$ is the roundtrip number such
that $E_{m}^{i}$ is the electric field of the $m$-th laser in the $i$-th
roundtrip. The coupling coefficient $K_{mn}$ is normalized by the self
coupling such that it depends only on the system geometry[25]. The effective
coupling strength is then determined by the simple equation:

$K(\theta)=K_{max}\sin[2](2\theta),$ (2)

where $K_{max}\approx-0.45$ is the numerically calculated coupling value for a
fully coupled cavity ($\theta=45^{o}$). It is worth noting that while this
configuration changes the spectral content of the cavity, it has a negligible
effect on phase locking for the parameters of this experiment[26, 27, 28].
Lastly, to achieve the best phase locking with a given coupling strength, we
use an intra-cavity adaptive optics method which reduces the effects of
uncontrolled frequency detuning in our system and improves our initial
conditions (see Fig. S6 in Supplemental Material)[29].

In the experiment, we introduce random frequency detuning to each laser in the
array, which limits the ability of the lasers to phase lock with each other
[30]. In each experimental realization, we generate a normally distributed
random frequency detuning pattern with a given standard deviation,
$\Omega_{RMS}$. When generating correlated detuning, we convolve the normally
distributed pattern with a two dimensional gaussian and renormalize it, such
that the resulting detuning pattern has a standard deviation of $\Omega_{RMS}$
and a correlation length of $\xi$. We perform 50 different realizations for
each value of $\Omega_{RMS}$ and $\xi$, and capture the resulting near field
(NF) and far field (FF) intensity distributions for each realization. The FF
intensity distribution ($I_{FF}$) is directly proportional to the Fourier
transform of the coherence function of the electric field[31], and so we
evaluate the phase locking in the array by calculating the average FF inverse
participation ratio (IPR), defined as:

$\text{IPR}=\frac{\int\int\differential{x}\differential{y}I_{FF}^{2}}{(\int\int\differential{x}\differential{y}I_{FF})^{2}},$
(3)

The IPR is commonly used to measure the localization of distributions and is
highly correlated with the average phase locked cluster size (see Fig. S2 in
Supplemental Material)[32, 33, 29]. Finally, we note that applying frequency
detuning and reducing phase locking introduces losses to the cavity due to the
dissipative nature of our coupling. To ensure that we operate in the same
regime throughout the experiment, i.e. far above the lasing threshold, we
modify the pumping power such that the total NF intensity stays roughly
constant (within $\pm 5\%$ of its initial value). A comparison between the
cavity performance in low and high detuning values is shown in Fig. S1 in the
supplemental material.

Results -  Fig. 2 shows the normalized IPR as a function of the detuning
spread, $\Omega_{RMS}/\absolutevalue{K}$, for normally distributed detuning
patterns ($\xi=0$).

Figure 2: Experimentally measured FF IPR as a function of the effective
detuning spread, $\Omega_{RMS}/\absolutevalue{K}$ for different coupling
values. IPR values are normalized between 1 and 0 and the dotted line
corresponds to $IPR=0.5$. Insets show the average FF intensity distribution at
$\tau_{c}\Omega_{RMS}=0,0.62$ for $\absolutevalue{K}=0.25$. Top right:
Approximate number of phase locked lasers as a function of detuning spread.
Dashed lines correspond to best fit to $y=\frac{a}{x^{b}+c}$.

Our results show, as expected, that increasing detuning disorder leads to a
smooth decrease in IPR (and hence phase locking). Furthermore, repeating the
experiment for different coupling values yields similar behavior, with the IPR
dropping to half its initial value for $\Omega_{RMS}/\absolutevalue{K}\approx
0.81(7)$, confirming that the behavior of the system is mostly determined by
the ratio between the coupling and detuning. The inset of Fig. 2 shows $N$,
the average phase locked cluster size, as estimated from the FF intensity
distribution [29], as a function of $\Omega_{RMS}$. While the value of $N$
should be identical for different coupling values for $\Omega_{RMS}=0$, we
find that it is not the case and that it instead differs significantly. This
is due to uncontrolled detuning, originating from system misalignments and
other errors, such that
$\Omega_{\text{actual}}=\Omega_{RMS}+\Omega_{\text{uncontrolled}}$. When
$\Omega_{RMS}\gg\Omega_{\text{uncontrolled}}$, the effects of the uncontrolled
detuning can be neglected. Using our results, we can estimate that
$\tau_{c}\Omega_{uncontrolled}\ll 0.1$ (additional related results can be
found in Fig. S6 in the Supplemental Material).

In the regime where the laser intensities are identical, our system is well
described by the Kuramoto model. Previous studies of Kuramoto oscillators with
nearest neighbor coupling have shown that for any finite coupling strength and
frequency detuning disorder, the maximal number of phase locked oscillators is
bounded [6, 34, 35]. In particular, for the case of a one-dimensional chain of
$N$ oscillators with nearest neighbor coupling $K$, the necessary condition
for phase locking is:

$\max{\absolutevalue{X_{j}}}=\max{\absolutevalue{\sum_{j=1}^{N}\Omega_{j}-\frac{1}{N}(\sum_{i=1}^{N}\Omega_{i})}}\equiv
K_{c}\leq K,$ (4)

where $\Omega_{j}$ is the frequency detuning of the $j$-th oscillator, and
$K_{c}$ is the critical coupling, which is the minimal required coupling for
phase locking. Simply stated, the condition for phase locking is that the
maximal accumulated detuning along the $N$ oscillators is smaller than the
coupling between two neighbors. In the case where $\Omega_{j}$ is a normally
distributed random vector, this expression is equivalent to the maximal
displacement of a random walker and implies that
$\max{\absolutevalue{X_{j}}}\propto\Omega_{RMS}\sqrt{N}$, such that
$N\propto\frac{K^{2}}{\Omega_{RMS}^{2}}$. Fitting the data from Fig. 2 to
$y=\frac{a}{x^{b}+c}$ we find values of $b=2.3(1),2.0(1),1.3(2)$ for the
different coupling values, with $c=0.044(12)$.

We now turn to focus on a single coupling strength, $\absolutevalue{K}=0.25$,
and study the effects of correlated detuning disorder (results for additional
coupling values can be found in Fig. S3-S5 the Supplemental Material). The
results are presented in Figs. 3-5. Fig. 3 shows the average FF IPR as a
function of $\Omega_{RMS}$ and the detuning correlation length $\xi$. It is
evident that for all values of $\xi$, the IPR decreases as $\Omega_{RMS}$ is
increased, as expected. However, the effect of correlated detuning disorder
appears to be non-trivial.

Figure 3: Average measured FF IPR as a function of $\Omega_{RMS}$ and the
correlation length $\xi$. The black dashed line shows the $IPR=0.4$ line. The
white dotted line highlights $\xi=8$, shown in detail in Fig. 4.

We highlight this behavior in Figs. 4-5: Fig. 4 shows the IPR as a function of
$\Omega_{RMS}$ for uncorrelated noise ($\xi=0$) and correlated detuning with
$\xi=8$ (highlighted with a white dotted line in Fig. 3). For
$\tau_{c}\Omega_{RMS}<0.31$, the IPR is lower (and phase locking is worse) for
the correlated detuning disorder compared to the uncorrelated case, while for
larger detuning values, the opposite is true. In order to validate the
results, we compare them to a numerical simulation of our system and find
qualitatively similar results (Fig. 4, right).

Figure 4: Experimentally measured (left) and numerically simulated (right) FF
IPR as a function of $\Omega_{RMS}$ for uncorrelated noise (black) and for
correlated noise with $\xi=8$ (green). Inset: Representative realizations of
noise, normalized to have the same RMS and shown on the same scale:
Uncorrelated (top) and with a correlation length $\xi=8$.

Fig. 5 shows the same data but as a function of $\xi$ for different values of
$\tau_{c}\Omega_{RMS}$. The figure reveals several phenomena: First, while for
low detuning values the highest IPR is measured at $\xi=0$, for high detuning
values the IPR monotonically increases with $\xi$. Second, the lowest IPR is
measured at an intermediate value, which changes as the detuning is increased.
Furthermore, the value of $\xi$ for which
$\text{IPR}(\Omega_{RMS},\xi)=\text{IPR}(\Omega_{RMS},0)$ decreases as we
increase $\Omega_{RMS}$. These effects are qualitatively replicated in
numerical simulations of our experiment.

Figure 5: Experimentally measured (left) and numerically simulated (right) FF
IPR as a function of correlation length $\xi$ for different values of
detuning, $\tau_{c}\Omega_{RMS}=0.12,0.23,0.36,0.62$ (blue, orange, yellow and
purple respectively). Dashed lines show the IPR value for $\xi=0$ for
reference.

Theoretical analysis and discussion -  In order to better understand our
results, we consider the locally coupled Kuramoto model in the case of
correlated disorder. The left panel of Fig. 6 shows the numerical evaluation
of Eq.(4) when assuming detuning disorder with gaussian correlations. We note
that for small cluster sizes, the critical coupling strength, which is equal
to the value of $\max{\absolutevalue{X_{j}}}$, is smaller when the correlation
length is longer - meaning it is harder to phase lock the same amount of
oscillators when the detuning disorder is uncorrelated compared to the
correlated case. However, for larger cluster sizes, the trend is reversed, and
the critical coupling grows larger with the correlation length of the
disorder. Notably, the crossing point between the two trends is roughly at
$N\sim\xi$. We derive an approximate analytical expression,
$\expectationvalue{\absolutevalue{X_{j}}}$:

$\displaystyle\expectationvalue{\absolutevalue{X_{j}}}^{2}$
$\displaystyle\approx\Omega_{RMS}^{2}\frac{\pi\xi
L}{8}\frac{1}{\frac{\xi}{\sqrt{2\pi}}(e^{-\frac{2L^{2}}{\xi^{2}}}-1)+L\operatorname{erf}(\frac{\sqrt{2}L}{\xi})}\times$
(5) $\displaystyle\sum_{i=1}^{L}$
$\displaystyle[\operatorname{erf}(\frac{j-i}{\xi})+\frac{1}{2}e^{-\frac{(j-i)^{2}}{\xi^{2}}}-\frac{j}{N}(\operatorname{erf}(\frac{N-i}{\xi})$
$\displaystyle+\frac{1}{2}e^{-\frac{(N-i)^{2}}{\xi^{2}}})+(1-\frac{j}{N})(\operatorname{erf}(\frac{i}{\xi})-\frac{1}{2}e^{-\frac{i^{2}}{\xi^{2}}})]^{2}$

The derivation of Eq.(5) can be found in the Supplemental Material. We compare
this analytical result and a numerical evaluation of Eq.(4) for the case of
correlated disorder in Fig. 6. It can be seen in the figure that the
analytical approximation captures the non-trivial dependence in the
correlation length of the disorder, namely that
$\expectationvalue{\absolutevalue{X_{j}}}$ decreases with $\xi$ for low N
values but increases with $\xi$ for higher $N$ values. analyzing the limiting
behavior of Eq.(5), reveals two distinct regimes: When $N\gg\xi$, i.e. the
phase locked cluster size is much larger than the correlation length of the
detuning disorder,
$\expectationvalue{\absolutevalue{X_{j}}}\rightarrow\sqrt{\xi N}$. This result
is equivalent to the displacement of a random walker which moves a distance of
$\xi>1$ for each step. In the opposite regime, $N\ll\xi$, i.e. the phase
locked cluster is much smaller than the correlation length of the disorder. In
this regime, $\expectationvalue{\absolutevalue{X_{j}}}\rightarrow
N\frac{N}{\xi}$, which is the result for evaluating Eq.(4) for a detuning
pattern of $\Omega_{i}=\sin(\frac{i}{\xi})$, meaning the disorder is
effectively approximated as a low frequency perturbation. A linear fit to the
log-log data of Fig. 6 for both regimes is shown in the bottom panel of the
figure. We find the best fit for $y=ax^{b}$ to be $b=-1.04(6),0.47(4)$ for the
$\xi\gg N,\xi\ll N$ regimes respectively, which is in good agreement with the
theoretical limiting behavior.

Figure 6: Left: Numerical evaluation of the Kuramoto model phase locking
condition for random detuning vectors with different correlation lengths. The
inset shows a zoomed-in view for low N values. Right: The analytical
approximation of the phase locking condition. Bottom: log-log plot of the
numerical data, shown as a function of $\xi$ for high (red) and low (blue)
values of $N$, the number of lasers in the phase locked cluster. Dashed lines
show the best fit to the colored-in points, which are in good agreement with
the analytical approximation.

The theoretical analysis provides us with a qualitative explanation of our
experimental results: The phase locked cluster size for the uncorrelated
detuning disorder, $N(\Omega_{RMS},K,\xi=0)$, determines a scale, $\xi_{c}$.
Since the phase locked cluster size is determined by the accumulated detuning
across the entire cluster, detuning disorder with correlation length
$\xi<\xi_{c}$ can be effectively considered as a correlated random walk, thus
accumulating a larger amount of detuning and causing a faster decay in phase
locking. In contrast, detuning disorder with $\xi>\xi_{c}$ is effectively a
low frequency perturbation in detuning along the cluster and causes a smaller
decay in phase locking. As we increase $\Omega_{RMS}$, $\xi_{c}$ becomes
smaller, until the maximal applied value of $\Omega_{RMS}=0.62$, where
$xi_{c}~{}0$ such that the minimal phase locking is measured for uncorrelated
disorder. This explanation is supported by both the numerical simulations of
our experiment and the theoretical analysis presented above. It is worth
noting that the theory on which we base our analysis is strictly one-
dimensional, in contrast to our two-dimensional system. Hence, a full
quantitative comparison between theory and experiment requires expanding the
current theoretical model to a higher dimension, which is outside the scope of
this work.

Conclusions -  In this work we experimentally measured the effects of
frequency detuning disorder on coupled laser arrays. Our results demonstrate
how increasing the amount of detuning disorder causes a gradual decay in phase
locking, which we have shown to depend on the ratio between the coupling
strength and the width of the detuning distribution. In addition, we have
shown that introducing a correlation length to the disorder yields a non-
monotonic behaviour: Correlated detuning disorder can be more or less damaging
to the phase locking of the array compared to uncorrelated disorder, depending
on the current number of phase locked lasers. We demonstrate this effect both
experimentally and numerically, and we show that they are in agreement with a
theoretical analysis based on the locally coupled Kuramoto model.

Our results provide an insight into the effects of disorder on many body
systems and can be expanded on in the future to different types of systems:
For example, long-range coupling can be expected to display significantly
different behavior, less dependent on the spatial distribution of the
disorder, where in the limit of all-to-all coupling correlation length can be
expected to have no effect on the system. In addition, by applying disorder in
a controlled way, as demonstrated in this work, it could be possible to
quantify the degree to which different systems are susceptible or resistant to
disorder, such as topologically protected states or intentionally introduced
defects.

Supplemental Materials

##  I Comparison of near field intensity distribution

As mentioned in the main text, the coupled laser array can be described by the
Kuramoto model equations when the lasers’ intensity is equal. However, the
introduction of frequency detuning to the lasers also induces increased losses
which are detuning dependent, caused both by the decrease in phase locking and
due to partial coupling between loss and detuning in our system (due to the
operation of the SLM). In order to avoid non-uniform lasing and to stay high
above the lasing threshold of the entire array, we increase our pump power
gradually to maintain roughly equal NF intensity. For each value of
$\Omega_{RMS}$ we measure the average total NF intensity, $I_{NF}^{tot}$, and
compare it to its value at $\Omega_{RMS}=0$. We increase or decrease the pump
power in order to keep $I_{NF}^{tot}$ within $5\%$ of its initial value. Fig.
S1 shows the NF and FF intensity distribution at the minimal and largest
$\tau_{c}\Omega{RMS}$ measured in the experiment. It can be seen that the NF
intensity is largely unaffected by the increased detuning, while FF shows
vastly reduced phase locking, evident in the wide diffraction peaks and strong
background.

Fig. S1: (a),(b): Normalized average NF intensity distributions for
$\tau_{c}\Omega_{RMS}=0,1.2rad$ respectively. (c),(d): Representative FF
intensity distributions.

##  II Quanifying phase locking with IPR

In order to validate our use of the IPR as a method for quantifying phase
locking, we perform numerical simulations of a coupled lasers array of $10x10$
lasers with varying random frequency detuning. The random detuning causes the
array splits into phase locked clusters, which decrease in size as the
detuning spread is increased [34, 35]. Fig. S2 shows a comparison between the
FF IPR and the average phase locked cluster size of the array. Our results
show that the two results are very well correlated, so it is reasonable to use
the FF IPR to quantify the average number of phase locked lasers.

Fig. S2: Left: IPR (blue) and normalized average phase locked cluster size
(red) calculated from LRE simulation, as a function of $\Omega_{RMS}$. Right:
IPR vs. phase locked cluster size. The data is highly correlated as shown by
the linear fit and correlation coefficient $r=0.99$.

##  III Correlated disorder results for weak coupling

We repeated the measurements shown in Fig. 3-5 in the main text with a weaker
coupling strength, $\absolutevalue{K}=\approx 0.12$. The results of both
experiments are shown side by side in Figures S3-S5. It is interesting to note
that the non-monotonic behavior in $\xi$ is absent from the weak coupling
data, and instead, increasing the correlation length of the disorder causes a
slower decay of the IPR. This is in agreement with our theoretical and
numerical analysis: Since the initial phase locked cluster size is smaller in
the case of weak coupling, the $\xi\ll\xi_{c}$ regime is unavailable. As a
result, all correlation lengths are of the scale of the cluster or larger,
hence cause a slower decay of the IPR compared to the uncorrelated case.

Fig. S3: Measured IPR as a function of $\Omega_{RMS},\xi$, for $|K|\approx
0.25\text{(left)},0.12\text{(right)}$ Fig. S4: Experimentally measured FF IPR
as a function of applied detuning, $\tau_{c}\Omega_{RMS}$ for
$\absolutevalue{K}=0.25$ (left) and $|K|=0.12$ (right). Different colored
plots correspond to different applied correlation length $\xi$. Fig. S5:
Experimentally measured FF IPR as a function of detuning disorder correlation
length $\xi$ for $\absolutevalue{K}=0.25$ (left) and $|K|=0.12$ (right).
Different colored plots correspond to different magnitudes of applied
detuning, $\Omega_{RMS}$.

##  IV The effect of intra-cavity adaptive optics

In previous work, we have developed an intra-cavity adaptive optics method
(AO) and shown its beneficial effect on phase locking[29]. In this experiment,
we use the same method in order to reduce the amount of uncontrolled frequency
detuning in the cavity. Fig. S6 shows the IPR as a function of detuning with
and without the application of adaptive optics. The figure shows that we
improve the initial phase locking our cavity by approximately $25\%$ as
measured by the IPR. Furthermore, it can be seen that for large
$\tau_{c}\Omega_{RMS}$ values, the two datasets coincide. This is likely due
to the fact that for these values,
$\Omega_{\text{uncontrolled}}\ll\Omega_{RMS}$, and hence the phase locking is
determined only by the change in $\Omega_{RMS}$. We note that this happens at
a relatively large value of $\tau_{c}\Omega_{RMS}$, hence highlighting the
importance of using our AO method to maximize cavity performance.

Fig. S6: Experimental measured FF IPR with (black) and without (red) adaptive
optics correction. Insets show the average FF intensity distribution at the
corresponding $\Omega_{RMS}$ value.

##  V Derivation of $\expectationvalue{\absolutevalue{X_{j}}}$

We first note some useful properties of normal random variables:

  1. 1.

Addition and Subtraction of random variables:

For two independent normally distributed random variables,

$Y=\mathcal{N}(\mu_{1},\sigma_{1}^{2});\;\;Z=\mathcal{N}(\mu_{2},\sigma_{2}^{2})$

Then the new random variable

$X=Y\pm Z$

is a normally distributed random variable,

$X=\mathcal{N}(\mu_{1}\pm\mu_{2},\sigma_{1}^{2}+\sigma_{2}^{2})$

  2. 2.

Correlation of normal variables: The correlation of two random variables is
defined by

$\rho(X,Y)=\frac{E[XY]-E[X]E[Y]}{\sigma_{X}\sigma_{Y}}=\frac{Cov(X,Y)}{\sqrt{Var[X]Var[Y]}}$

Where

$\sigma_{x}=\sqrt{E[X^{2}]-E[X]^{2}}$

And results in $-1\leq\rho\leq 1$

  3. 3.

Addition of correlated variables:

In general, for any set of random variables,

$Var[\sum
X_{i}]=\sum_{i}\sum_{j}Cov(X_{i},X_{j})=\sum_{i}Var[X_{i}]+2\sum_{i}\sum_{j\geq
i}Cov(X_{i},X_{j})$

For normally distributed variables then,

$X=Y\pm
Z=\mathcal{N}(\mu_{1}\pm\mu_{2},\sigma_{Y}^{2}+\sigma_{Z}^{2}\pm\rho\sigma_{Y}\sigma_{Z})$

By extension, we see that when we take a multiple of a single variable,

$Var[aX_{i}]=aVar[X_{i}]+2\matrixquantity(a\\\
2)Cov[X_{i},X_{i}]=a\sigma_{X}^{2}+a(a-1)\sigma_{X}^{2}=a^{2}\sigma_{X}^{2}$

so if

$X=aY$ $X=\mathcal{N}(a\mu_{Y},a^{2}\sigma_{Y}^{2})$

Lastly, we note that

$Cov[aX,bY]=E[abXY]=abE[XY]$

We turn to derive an expression for the expression of Eq. (5):

$\expectationvalue{\absolutevalue{X_{j}(N,\xi)}}=\expectationvalue{\absolutevalue{\sum_{i=1}^{j}(\Omega_{i}-\bar{\Omega})}}$

We consider a series of $L$ independent normal random variables,

$Z_{i}=\mathcal{N}(0,\sigma_{Z}^{2})$

And construct a series of correlated random variables by convolving the $Z$
variables with a gaussian with waist $\xi$:

$r_{j}=\sum_{i=1}^{L}e^{-\frac{(i-j)^{2}}{\xi^{2}}}Z_{i}$

We note that $X_{j}$ can be expressed in terms of partial sums of the series,

$X_{j}(N,\xi)=\sum_{i=1}^{j}(Z_{i}-\bar{Z})=S_{j}(\xi)-\frac{j}{N}S_{N}(\xi)$

And so we start by calculating their value:

$S_{n}=\sum_{j=1}^{n}r_{j}=\sum_{j=1}^{n}\sum_{i=1}^{L}e^{-\frac{(i-j)^{2}}{\xi^{2}}}Z_{i}$
$X_{n}=S_{n}-\frac{n}{N}S_{N}=\sum_{i=1}^{L}[\sum_{j=1}^{n}e^{-\frac{(i-j)^{2}}{\xi^{2}}}Z_{i}-\frac{n}{N}\sum_{k=1}^{N}e^{-\frac{(i-k)^{2}}{\xi^{2}}}Z_{i}]$

we approximate the series as an integral using the first order Euler-Maclaurin
approximation:

$\approx\sum_{i=1}^{L}[\int_{0}^{n}e^{-\frac{(j-i)^{2}}{\xi^{2}}}\differential{j}+\frac{e^{-\frac{(n-i)^{2}}{\xi^{2}}}-e^{-\frac{i^{2}}{\xi^{2}}}}{2}-\frac{n}{N}\int_{0}^{N}e^{-\frac{(i-k)^{2}}{\xi^{2}}}\differential{k}-\frac{n}{N}\frac{e^{-\frac{(N-i)^{2}}{\xi^{2}}}-e^{-\frac{i^{2}}{\xi^{2}}}}{2}]Z_{i}$
$=\sum_{i=1}^{L}[\operatorname{erf}(\frac{n-i}{\xi})+\operatorname{erf}(\frac{i}{\xi})+\frac{e^{-\frac{(n-i)^{2}}{\xi^{2}}}-e^{-\frac{i^{2}}{\xi^{2}}}}{2}-\frac{n}{N}\operatorname{erf}(\frac{N-i}{\xi})-\frac{n}{N}\operatorname{erf}(\frac{i}{\xi})-\frac{n}{N}\frac{e^{-\frac{(N-i)^{2}}{\xi^{2}}}-e^{-\frac{i^{2}}{\xi^{2}}}}{2}]\frac{\sqrt{\pi}}{2}\xi
Z_{i}$

While additional orders may refine the result, we claim that this first order
approximation is sufficient for the purpose of this work.

For $n\rightarrow N$, we expect $X_{n}\rightarrow 0$, and indeed we see that
it is the case.

After creating a series of correlated variables, we normalize it to have a
unity standard deviation. For that reason, we need to calculate the standard
deviation of $L$ correlated variables. For a single variable, we get:

$E[r_{n}^{2}]=E[(\sum_{i=1}^{L}e^{-\frac{(i-n)^{2}}{\xi^{2}}}Z_{i})^{2}]=\sum_{i=1}^{L}e^{-\frac{2(i-n)^{2}}{\xi^{2}}}E[Z_{i}^{2}]=\frac{\sqrt{\pi}}{2\sqrt{2}}\xi\sigma_{Z}^{2}[\operatorname{erf}(\frac{\sqrt{2}n}{\xi})+\operatorname{erf}(\frac{\sqrt{2}(L-n)}{\xi}]$

where we used the fact that the $Z_{i}$ variables are independent, so
$E[Z_{i}Z_{j}]=\delta_{ij}$. And so for the entire series of $L$ variables, we
get:

$\sigma_{L}^{2}=\frac{\sqrt{\pi}}{2\sqrt{2}}\frac{\xi\sigma_{Z}^{2}}{L}\sum_{n=1}^{L}[\operatorname{erf}(\frac{\sqrt{2}n}{\xi})+\operatorname{erf}(\frac{\sqrt{2}(L-n)}{\xi}]\approx\frac{2\xi\sigma_{Z}^{2}}{L}\int_{0}^{L}\operatorname{erf}(\frac{\sqrt{2}{n}}{\xi})\differential{n}$

We note here a useful integral solution,

$\int_{0}^{x}a\operatorname{erf}(\frac{(x^{\prime}-b)}{a})\differential{x}=\frac{a^{2}}{\sqrt{\pi}}(e^{-\frac{(x-b)^{2}}{a^{2}}}-e^{-\frac{b^{2}}{a^{2}}})+a(x-b)\operatorname{erf}(\frac{x-b}{a})-ab\operatorname{erf}(\frac{b}{a})$

such that

$\sigma_{L}^{2}=\frac{2\xi\sigma_{Z}^{2}}{L}[\frac{\xi}{\sqrt{2\pi}}(e^{-\frac{2L^{2}}{\xi^{2}}}-1)+L\operatorname{erf}(\frac{\sqrt{2}L}{\xi})]$

Finally, using the properties of normal variables shown above, we get that:

$\expectationvalue{\absolutevalue{X_{n}}}^{2}\approx\frac{\pi\xi
L}{8}\frac{\sum_{i=1}^{L}[\operatorname{erf}(\frac{n-i}{\xi})+\frac{1}{2}e^{-\frac{(n-i)^{2}}{\xi^{2}}}+(1-\frac{n}{N})(\operatorname{erf}(\frac{i}{\xi})-\frac{1}{2}e^{-\frac{i^{2}}{\xi^{2}}})-\frac{n}{N}(\operatorname{erf}(\frac{N-i}{\xi})]^{2}-\frac{1}{2}e^{-\frac{(N-i)^{2}}{\xi^{2}}})}{\frac{\xi}{\sqrt{2\pi}}(e^{-\frac{2L^{2}}{\xi^{2}}}-1)+L\operatorname{erf}(\frac{\sqrt{2}L}{\xi})}$

## References

  * [1] Thomas G Kiely and JK Freericks.  Relationship between the transverse-field ising model and the x y model via the rotating-wave approximation.  Physical Review A, 97(2):023611, 2018.
  * [2] Malte Henkel.  Statistical mechanics of the 2d quantum xy model in a transverse field.  Journal of Physics A: Mathematical and General, 17(14):L795, 1984\.
  * [3] Róbert Juhász, István A Kovács, and Ferenc Iglói.  Random transverse-field ising chain with long-range interactions.  Europhysics Letters, 107(4):47008, 2014.
  * [4] Ross H McKenzie.  Exact results for quantum phase transitions in random xy spin chains.  Physical review letters, 77(23):4804, 1996.
  * [5] Daniel S Fisher.  Random transverse field ising spin chains.  Physical review letters, 69(3):534, 1992.
  * [6] Juan A Acebrón, Luis L Bonilla, Conrad J Pérez Vicente, Félix Ritort, and Renato Spigler.  The kuramoto model: A simple paradigm for synchronization phenomena.  Reviews of modern physics, 77(1):137, 2005.
  * [7] Kurt Wiesenfeld, Pere Colet, and Steven H Strogatz.  Frequency locking in josephson arrays: Connection with the kuramoto model.  Physical Review E, 57(2):1563, 1998.
  * [8] Brad R Trees, Vinod Saranathan, and David Stroud.  Synchronization in disordered josephson junction arrays: Small-world connections and the kuramoto model.  Physical Review E, 71(1):016215, 2005.
  * [9] MK Stephen Yeung and Steven H Strogatz.  Time delay in the kuramoto model of coupled oscillators.  Physical review letters, 82(3):648, 1999.
  * [10] Naotomo Takemura, Kenta Takata, Masato Takiguchi, and Masaya Notomi.  Emulating the local kuramoto model with an injection-locked photonic crystal laser array.  Scientific Reports, 11(1):8587, 2021.
  * [11] Shir Shahal, Ateret Wurzberg, Inbar Sibony, Hamootal Duadi, Elad Shniderman, Daniel Weymouth, Nir Davidson, and Moti Fridman.  Synchronization of complex human networks.  Nature communications, 11(1):3854, 2020.
  * [12] Vishwa Pal, Chene Tradonsky, Ronen Chriki, Asher A. Friesem, and Nir Davidson.  Observing dissipative topological defects with coupled lasers.  Phys. Rev. Lett., 119:013902, Jul 2017.
  * [13] Lechen Yang, Guangrui Li, Xiaomei Gao, and Ling Lu.  Topological-cavity surface-emitting laser.  Nature Photonics, 16(4):279–283, 2022.
  * [14] Alex Dikopoltsev, Tristan H. Harder, Eran Lustig, Oleg A. Egorov, Johannes Beierlein, Adriana Wolf, Yaakov Lumer, Monika Emmerling, Christian Schneider, Sven Höfling, Mordechai Segev, and Sebastian Klembt.  Topological insulator vertical-cavity laser array.  Science, 373(6562):1514–1517, 2021.
  * [15] Christian E Rüter, Konstantinos G Makris, Ramy El-Ganainy, Demetrios N Christodoulides, Mordechai Segev, and Detlef Kip.  Observation of parity–time symmetry in optics.  Nature physics, 6(3):192–195, 2010.
  * [16] Geva Arwas, Sagie Gadasi, Igor Gershenzon, Asher Friesem, Nir Davidson, and Oren Raz.  Anyonic-parity-time symmetry in complex-coupled lasers.  Science advances, 8(22):eabm7454, 2022.
  * [17] Micha Nixon, Eitan Ronen, Asher A. Friesem, and Nir Davidson.  Observing geometric frustration with thousands of coupled lasers.  Phys. Rev. Lett., 110:184102, May 2013.
  * [18] Chene Tradonsky, Simon Mahler, Gaodi Cai, Vishwa Pal, Ronen Chriki, Asher A. Friesem, and Nir Davidson.  High-resolution digital spatial control of a highly multimode laser.  Optica, 8(6):880–884, Jun 2021.
  * [19] Hui Cao, Ronen Chriki, Stefan Bittner, Asher A Friesem, and Nir Davidson.  Complex lasers with controllable coherence.  Nature Reviews Physics, 1(2):156–168, 2019.
  * [20] J. A. Arnaud.  Degenerate optical cavities.  Appl. Opt., 8(1):189–196, Jan 1969.
  * [21] Chene Tradonsky, Vishwa Pal, Ronen Chriki, Nir Davidson, and Asher A Friesem.  Talbot diffraction and fourier filtering for phase locking an array of lasers.  Applied Optics, 56(1):A126–A132, 2017.
  * [22] James R. Leger.  Lateral mode control of an algaas laser array in a talbot cavity.  Applied Physics Letters, 55(4):334–336, 1989.
  * [23] Eitan Ronen and Amiel A. Ishaaya.  Phase clusters induced by degeneracy in a phase locked fiber laser array.  IEEE Journal of Quantum Electronics, 47(12):1526–1530, 2011.
  * [24] David Mehuys, William Streifer, Robert G. Waarts, and David F. Welch.  Modal analysis of linear talbot-cavity semiconductor lasers.  Opt. Lett., 16(11):823–825, Jun 1991.
  * [25] Fabien Rogister, K. Scott Thornburg, L. Fabiny, Michael Möller, and Rajarshi Roy.  Power-law spatial correlations in arrays of locally coupled lasers.  Phys. Rev. Lett., 92:093905, 2004.
  * [26] Sharona Sedghani Cohen, Vardit Eckhouse, Asher A Friesem, and Nir Davidson.  Single frequency lasing using coherent combining.  Optics communications, 282(9):1861–1866, 2009.
  * [27] Moti Fridman, Micha Nixon, Eitan Ronen, Asher A Friesem, and Nir Davidson.  Phase locking of two coupled lasers with many longitudinal modes.  Optics letters, 35(4):526–528, 2010.
  * [28] P Smith.  Stabilized, single-frequency output from a long laser cavity.  IEEE Journal of Quantum Electronics, 1(8):343–348, 1965.
  * [29] Amit Pando, Sagie Gadasi, Asher Friesem, and Nir Davidson.  Improved laser phase locking with intra-cavity adaptive optics.  Opt. Express, 31(4):6947–6955, Feb 2023.
  * [30] S. Yu. Kourtchatov, V. V. Likhanskii, A. P. Napartovich, F. T. Arecchi, and A. Lapucci.  Theory of phase locking of globally coupled laser arrays.  Phys. Rev. A, 52:4089–4094, Nov 1995.
  * [31] Leonard Mandel and Emil Wolf.  Coherence properties of optical fields.  Reviews of modern physics, 37(2):231, 1965.
  * [32] Dirk Witthaut and Marc Timme.  Kuramoto dynamics in hamiltonian systems.  Phys. Rev. E, 90:032917, Sep 2014.
  * [33] Stefano Longhi.  Non-hermitian laser arrays with tunable phase locking.  Optics Letters, 47(8):2040–2043, 2022.
  * [34] Steven H Strogatz and Renato E Mirollo.  Collective synchronisation in lattices of nonlinear oscillators with randomness.  Journal of Physics A: Mathematical and General, 21(13):L699, 1988\.
  * [35] Hidetsugu Sakaguchi, Shigeru Shinomoto, and Yoshiki Kuramoto.  Local and grobal self-entrainments in oscillator lattices.  Progress of Theoretical Physics, 77(5):1005–1010, 1987.