Emergent charge order near the doping-induced Mott-insulating quantum phase transition in Sr$_{3}$Ru$_{2}$O$_{7}$

We investigate the double layered Sr$_{3}$(Ru$_{1-x}$Mn$_{x}$)$_{2}$O$_{7}$ and its doping-induced quantum phase transition (QPT) from a metal to an antiferromagnetic (AFM) Mott insulator. Using spectroscopic imaging with the scanning tunneling microscope (STM), we visualize the evolution of the electronic states in real- and momentum-space. We find a partial-gap in the tunneling density of states at the Fermi energy (E$_{F}$) that develops with doping to form a weak Mott insulating ({\Delta} ~ 100meV) state. Near the QPT, we discover a spatial electronic reorganization into a commensurate checkerboard charge order. These findings share some resemblance to the well-established universal charge order in the pseudogap phase of cuprates. Our experiments therefore demonstrate the ubiquity of the incipient charge order that emanates from doped Mott insulators.

symmetry states emerge near the QPT. Figure 1(c) shows a topographic STM image of single crystal Sr 3 (Ru 1-x Mn x ) 2 O 7 with x=5% cleaved in-situ in our ultrahigh-vacuum STM exposing the SrO surface. The Fourier transform (FT) of the topograph reveals atomic Bragg peaks at the corners of the map corresponding to (±2/a, 0) and (0, ±2/a) ( Fig.1(d)). Note that, throughout this paper, we use tetragonal notation with lattice constant a ~ 3.9Å. Additional satellite reflections at (±/a, ±/a) correspond to the doubling of the bulk unit cell owing to the rotational distortion of the RuO 6 octahedra ( Fig.1(a),(d)) as previously observed [11,20]. The Mn dopants can be identified by their bias-dependent signatures in the topographic images ( Fig.1(c)). They can be best visualized at a sample bias V ~ 1 V as bright spots for different doped samples ( Fig.1(e)). Counting the number of Mn dopants on the surface reveals a systematic down shift (~ 2%) with a lower surface doping as compared to the expected bulk (see Supplementary Information).
The local impact of Mn-dopants on the spatial electronic structure, over a range of doping spanning the QPT, can be visualized through the spectroscopic imaging shown in Fig.2. At low Mn concentration (x=1%; much below the QPT) the local electronic density of states (LDOS) reveal a particle-hole asymmetric V-shaped partial-gap of ~ ±30 meV near the Fermi energy at T = 13 K ( Fig.2(a)-(c)). The density of states is similar to what has been previously observed in the undoped parent compound [24]. In particular, a small hump at E F previously identified as a van-Hove singularity [9,24], is also present in our spectra obtained away from Mn dopants (Fig.2(c) and supplementary information).
The effect of the Mn, at this low concentration, shows modifications of the local gap predominantly near E F (see Fig.2(b),(c) and see supplementary information). Conductance maps dI/dV (r,V) ≡ g(r, V) locally at locations r and at a constant sample bias V, displayed in Fig.2(b), reveal that, in this dilute limit the Mn-dopants spatially modify the LDOS within a length scale of 20Å.
Increasing the Mn-doping towards the QPT reveals a striking spatial reorganization of the local electronic states. Figure 2(d)-(f) display a topographic image and the LDOS, showing the evolution of the partial-gap in the x=5% sample. Conductance maps g(r, V) reveal a bi-directional long-range checkerboard charge modulations (Fig.2(e)); particularly V=-40meV and V=-100meV maps) embedded within the long-wavelength inhomogeneity. The absence of this charge ordering in the real space images of the x=1% sample ( Fig.2(b)), as well as in earlier STM experiments on the undoped [11,24] and Mndoped compound at 100 K [20], indicate its origin not to be related to any surface reconstruction, rather to an electronic instability induced by Mn-doping at low temperatures. Further indication of the electronic nature of the ordering phenomenon can be concluded by the absence of any structural peaks, corresponding to the observed modulations, in neutron scattering [14], in recent low energy electron diffraction (LEED) [20], and in non-resonant hard x-ray scattering [18]. The observed charge modulations, along the Ru-O bond direction with a periodicity of 2a (corresponding to a wave vector q = (±/a, 0) and (0, ±/a) ≡ Q*) therefore corroborate an emergent electronic charge ordering near the metalto-Mott-insulator phase transition.
At even higher doping (x=10%), the STM spectra reveals strong inhomogeneity with areas of Vshaped partial gaps mixed with a complete, yet inhomogeneous, gapping of the Fermi surface into a Ushaped insulating gap of ~100meV ( Fig.2(i) and supplementary). This is indicative of a close proximity to the metal-to-Mott insulator transition [25]. The evolution of the gap with doping resembles the opposite trend of what is seen in the hole-doped, spin-orbit driven, Mott insulator Sr 3 (Ir 1-x Ru x ) 2 O 7 [25]. The conductance maps g(r, V) ( fig.2(h)) reveal inhomogeneous short-range modulations, equally populated along the two atomic directions with their periodicity locked to 2a. The glassy nature of the observed charge order could be a consequence of the larger electronic inhomogeneity introduced by the Mndopants and the transition into the Mott insulating state.
It is important to emphasize that the doping-induced metal to insulator transition evaluated from the surface LDOS at T = 13 K (Blue points in Fig.1(b)) closely follows that extracted from bulk resistivity measurements (gray area in Fig.1(b)) [16]. The slight offset may be a result of the lower doping on the surface as compared to the bulk. This may also explain the slight enhancement of the octahedral distortion on the surface of Sr 3 (Ru 1-x Mn x ) 2 O 7 as compared to the bulk [26].
To visualize the momentum-space structure of these modulations and their doping dependence, we display in Fig.3 the real-space conductance maps and selected FTs at several energies taken on a larger surface area for optimal momentum-space resolution. For x=1%, in addition to the Bragg peaks at (±2/a, 0) and (0, ±2/ a) and the octahedral rotation satellite peaks at (±/a, ±/a), the FTs show dispersive features (q 2 , q 7 and q 8 ) along the high symmetry directions, corresponding to quasiparticle interference ( Fig.3(a), (b)) (here we follow the same notations as reference [11], similar to those previously seen in the STM studies of 1% Ti-doped Sr 3 Ru 2 O 7 [11]. The extracted dispersions along (±/a, ±/a) are shown in Fig.3(g). While no clear peak is observed at the Q* location in this 1% sample, there is a faint and broad intensity in the FT around this wave vector (q 1 in Fig.3b, whose dispersion is difficult to extract due to its proximity to q7 and their weak intensity), which locally originates only from areas with Mn dopants (more detail in supplementary information).
In contrast, for x=5%, the FTs indicate distinct sharp and non-dispersive peaks at Q* = (± /a, 0) and (0, ± /a) that correspond to a real-space modulation with a periodicity of 2a ( Fig.3(d),(h)). In addition to the Q* modulations, we also observe strong non-dispersive modulations coinciding with the AFM ordering wavevector, Q AFM = (±/2a, ±/2a) (Fig.3(d),(g)). In principle, without a spin-polarized tip, STM is not sensitive to a magnetic contrast. Our experiments indicate that the AFM order also has a charge component, which was suggested by previous REXS experiments [16]. Moving into the Mott state (x=10%) reveals the Q* modulations to broaden, yet dominate and remain locked to the same commensurate momentum, independent of doping ( Fig.3(f),(h)).
These experiments, unsurprisingly, indicate strong interplay between spin and charge instabilities.
Close to the QPT (x ~ 5%) the system is susceptible to long-range charge ordering at Q* and Q AFM (Fig.3) as well as a competing short-range AFM ordering at Q AFM (as revealed by previous REXS [16,18]). Signs of such a competition can also be seen in the suppressed AFM correlation length with decreasing temperature in REXS [16,18]. This presents an exotic electronic phase in which charge and spin orders simultaneously coincide at the same wavevector. With further doping, however, the charge ordering dominates at Q*, whereas a long-range E-type AFM order (with negligible charge contribution, as seen in our STM data ( Fig.3(f)) dominates at Q AFM [14,18]. It is important to note that although STM is primarily sensitive to charge order, some forms of orbital orders can also be detected through the STM tip sensitivity to different orbital wave functions [27,28] and therefore orbital contribution to Q AFM cannot be omitted as a possible scenario.
To provide insights into the origins of the emergent broken symmetry states at Q* and Q AFM , we compute the bare susceptibility (q,) using a 12 orbital tight-binding model [29], fit to match the Fermi surface measured by ARPES on Sr 3 Ru 2 O 7 (see methods and supplementary information). The calculations of the bare susceptibility, (q,0), whose peaks correspond to the wavevectors with strong Fermi surface nesting, reveal broad features around Q* = (/a, 0) and Q AFM = (/2a, /2a) , which are enhanced with increased Mn-doping (see supplementary information). As previously shown, the effect of Mn-dopants on the Sr 3 Ru 2 O 7 band structure can be partly accounted for by the suppression of the RuO 6 octahedral rotation in the tight-binding band structure [31]. As doping increases, the suppressed rotational-distortion (modeled as suppression in the hopping parameter [31]) reduces the coupling between the d xz and d yz orbitals making the  2 band more square-like. This itself enhances the nesting at Q * and Q AFM , leading to an increase in (q,0) (see supplementary information). Our theoretical analysis therefore indicates that a broad Fermi surface nesting does exist near Q* and Q AFM near x ~ 5%. However, we would like to emphasize that nesting alone cannot explain the emergent charge order at Q* and it will require enhanced anisotropic electron-phonon (e-p) coupling and/or strong coulomb correlations to induce a commensurate charge ordering instability. Recently, strong e-p coupling have been observed at the charge ordering wavevectors in YBa 2 Cu 3 O 6.6 [32], Bi 2 Sr 2 CaCu 2 O 8+x [33] and NbSe 2 [34]. In addition, a significant softening in the temperature dependence of the B 1g phonon mode has also been observed to coincide with the metal to Mott insulator transition in Ca 3 Ru 2 O 7 pointing towards the possibility of similar physics in Mn-doped Sr 3 Ru 2 O 7 . Inelastic resonant x-ray or neutron scattering [35,36] can probe the momentum dependence of the e-p coupling and such future experiments may offer a more complete picture of the origin of the observed instabilities in Sr 3 (Ru 1-x Mn x ) 2 O 7 .
Finally, we examine the structure and symmetry of the charge ordering at Q*. Figure 4(a) shows a zoom-in of the conductance map for the x=10% sample (white box in Fig.3(e)). While the charge order can be clearly observed, the long-wavelength inhomogeneity in the conductance map obscures the detail of the ordering patterns. To visualize the real space character of these patterns, we Gaussian-filter the Q* peaks in the raw unsymmetrized FTs (see supplementary information) and inverse Fourier transform back into real space. Figure 4(b) shows the charge order map corresponding to Fig.4(a). The figure displays a mixture of spatially inhomogeneous checkerboard and stripe patterns. To quantitatively distinguish between the two glassy orders, we compute the amplitudes of the two orthogonal components of the charge order ( Fig.4(c), (d)) and calculate their cross-correlation ( Fig.4(e)) (see Supplementary   Information). In a checkerboard order, the two orthogonal components are locally correlated (positive cross-correlation). A stripe order leads to a local anti-correlation between the two orthogonal charge ordering amplitudes. The positive cross-correlation observed as a function of bias ( Fig.4(e)) indicates globally a four-fold symmetric bi-directional charge order.
In conclusion, our experiments reveal non-dispersive commensurate charge order instability to emerge near the metal to AFM Mott insulator transition in Sr 3 Ru 2 O 7 . These findings draw similarities to the ubiquitous charge order seen in STM and x-ray experiments on cuprates [37][38][39][40][41][42][43][44]. While the metal to Mott insulator transition in Sr 3 Ru 2 O 7 , which is likely related to structural distortion and enhanced Coulomb correlations by 3d Mn replacing 4d Ru, is quite different than the charge doping-induced Mott insulator to metal transition in cuprates, it is quite intriguing that a charge ordered state with similar realspace patterns exists near the Mott transition of both systems. Finally, similar charge modulation patterns on strong spin-orbit coupled Sr 2 IrO 4 have recently been found near its doping-induced Mott insulator to metal transition with a similar periodicity to our Q* [45]. Overall, all these experiments, independent of their detailed intrinsic electronic structure, point to ubiquitous tendencies towards some form of charge order emerging at the border of Mott insulators.

Sample preparation and characterization:
The single crystals of Sr 3 (Ru 1-x Mn x ) 2 O 7 used in our experiments were grown at RIKEN laboratory using the floating zone technique. These single crystals were obtained under 1M Pa. atmosphere of argon gas containing 10% oxygen. Small flat samples were cut down to ~ (2 x 2 x 0.5mm 3 ) and attached to aluminum plates using H20E conducting epoxy. Aluminum cleaving posts were then attached to samples perpendicular to the a-b cleaving plane using H74F non conducting epoxy. Samples were cleaved at room temperature in ultrahigh vacuum by knocking the post off, and then immediately transferred in situ to the custom RHK Pan STM head, which had been cooled down to T ~ 13 K. STM topographs were taken in constant current mode and dI/dV measurements were performed using a standard lock-in technique with a reference frequency of 0.921 kHz. PtIr tips were used in all experiments. Tips were prepared prior to each experiment on a Cu(111) surface that had been treated with several rounds of sputtering and annealing, and then placed into the microscope head to cool. The sample and Cu were placed next to each other inside the microscope head so as to reduce exposure and preserve tip structure when moving from one to the other. The data presented in this paper were collected from approximately 3 successfully cleaved samples for each doping percentage (1%, 5%, 7.5%, and 10%). There were negligible differences measured between different cleaves and between different areas on a sample.

Tight binding calculations:
To calculate the bare susceptibility, we start from a 12 orbital tight-binding model that contains t 2g orbitals, the bilayer structure, and the spin-orbit coupling. The hopping parameters are determined by fitting the Fermi surfaces observed in ARPES. Details can be found in Ref. [29]. The resulting Hamiltonian in the momentum space can be expressed as where Φ = ( , , , , , , ,− , , + , , , + , , , + ,− , ) are the electron annihilation operators, k is the in-plane momentum, k z is the out-of-plane momentum associated with the bilayer structure, Q=(, ) is the wavevector associated with the octahedral rotation, and s is the spin index. The matrix element h k (k z ) can be found in Ref. [29]. The retarded bare susceptibility in a multiorbital system can in general be written as a tensor of Finally, we compute ( , 0) = ∑ 1 , 2 , 3 , 4 ( , 0) 1 , 2 , 3 , 4 which exhibits peaks at wavevectors corresponding to the Fermi surface nesting.
As the doping is introduced, the octahedral rotation is reduced. The leading effect of this octahedral rotation is to enable hopping between different orbitals on nearest neighbor sites. The details of this model is discussed in Ref. [29]. The reduced octahedral rotation results in less hybridization between the d xz and d yz orbitals. Note that the current model cannot capture the insulating behavior in the paramagnetic state. showing the doping-temperature (right-axis) transition from a PM metal to an AFM Mott insulator [16].
The blue circles represent the doping dependence of the STM LDOS at E F measured at T=13 K (left axis) normalized to value at 200meV. The blue square is extracted from ref. [24].

Identifying the Mn dopants and density:
Identifying exact Mn dopant sites is important to many aspects of our analysis. The sample cleaves along the Sr-O plane revealing a layer with the Sr atoms appearing as a bright lattice of atoms and the apical Oxygen and sublayer Ru appearing as the lower intensity areas between the Sr lattice. The Mn atoms replace the Ru, which means they are expected to appear in between the surface atoms (Fig. S1). Figure S1 shows the impurity sites appearing as a bright spot centered over the Ru/Mn depression that is in between the surface Sr atoms as expected.  There is a strong bias dependence on the appearance of the Mn dopants. Figure S2(a)-(e) shows several topographies taken on a 10% Mn doped sample, on the same field of view, with increasing bias. The Mn impurities appear as bright spots centered in between the surface atoms. As the bias is increased, the intensity of the bright spots increases until they dominate the image around 1V. Figure S2(f) shows a similar appearance for the Mn impurity at 1V in the 1% sample and a comparison between e and f shows how the impurity density increases between the 1% and 10% samples.
Counting these impurities on several topographies of varying sizes reveals a Mn percentage of (1.3 ± 0.1)% for the 1%, (2.6 ± 0.5)% for the 5%, (5.2 ± 0.4)% for the 7.5%, and (7.8 ± 0.4)% for the 10%. This shift in Mn % from the expected bulk (The bulk Mn composition in our samples was checked by inductively coupled plasma (ICP) analysis [1].) to the surface is shown in Fig. S3. Similar differences between the surface and the bulk were reported [2] for independently grown samples. The distinct size, appearance, number, and location of the impurity sites allows us to differentiate them from other defects intrinsic to the sample.

Determination of Gap Width:
FIG. S4: Determination of gap width. Shown in top panel are spectra for a 1%, a fully gapped 10% and a partial gapped 10% taken from the spectra used in fig. 2 of the main text. The bottom panel shows the 2 nd derivative of the spectra shown above. The size of the gap can be determined by the local minimum of the second derivative and are marked by the blue area for the fully gapped 10% spectra. All gap sizes used in this work are determined using this method.

Local effects of Mn impurities on x=1% sample:
The local effect of the Mn impurities on the spectra in the 1% doped sample is a local suppression of the density of states near EF in the vicinity of the van-Hove singularity on impurity sites. Figure S5 (a) shows a topography on a 1% sample with impurity locations marked. The purple circles are the Mn sites described above, and the red circles are defect sites. Spectra on and off impurity sites is shown in Fig S5 (b),(c), revealing a strong feature near EF on clean areas that is suppressed near impurities. Figure S5 (d), shows a spatial map of the intensity of the hump at EF (van-Hove map) obtained by plotting the conductance map corresponding to the energy (2*G(r, 2mV) -G(r,-6.5mv) -G (r, 8.5mV)). An overlay of impurity sites reveals the impurities to be grouped on areas of suppressed van-Hove peak, with a strong van-Hove signal only appearing on areas free of impurities. While there is no charge order and no strong Q* peak present in the 1% sample, a weak feature (q1) is still present at the same wavevector. Following a similar procedure described below for the Q* peak, Fig S5 (f), shows a spatial map of the amplitude of this broad feature, with dark areas corresponding to high intensity and light areas to low or zero intensity. Overlaying impurity locations reveals strong intensity locations (corresponding to broad modulations with the q1 wavevector) to correspond to areas of Mn impurities, with zero contribution coming from clean areas, indicating this feature to originate from local effects of the Mn impurities. Cross correlation taken between the van-Hove map and the q1 amplitude map gives a value of -0.43. This large negative value (anticorrelation) indicates that where the van-Hove peak is strong (away from Mn impurities), the q1 amplitude is weak. This analysis indicates that whether the origin of the q1 feature is related to quasiparticle interference off impurities or the beginning of local short-range ordering, it originates from the Mn dopants. In addition, no ordering can be observed locally on areas away from Mn establishing the absence of long range order.

Raw and Symmetrized FT of STM conductance maps
FIG. S6: Raw and Symmetrized FT of STM conductance maps. Raw Fourier Transform (FT) of the conductance maps for 5% and 10% Mn samples shown in Fig 3 of the main text (left). To enhance the noise to signal ratio, the FTs are Gaussian smoothed (middle) and four-fold symmetrized (right). The bright spot in the center (very small q) due to long wavelength fluctuations is removed by subtracting Gaussian-filtered conductance map from the original raw conductance value at every pixel. The 2D Fourier transform gets rid of only the low q-vector components. The symmetrized FTs are properly rotated to have the Bragg peaks at the corners of the image. Figure S7 shows line cuts taken along the (π/a, π/a) direction for the symmetrized 1% data shown in Fig 3, with their corresponding fits shown as a dashed red line. According to [3], two QPI features are expected to be observed here, namely q2 and q8. Our linecuts show two sets of dispersive features which we attribute to these two wave vectors. The q8 and q2 features were fit separately and then combined in the figure. The fits are constructed from a background made of a 2 nd order polynomial. The q8 feature is fit to a single Gaussian with the q value extracted from the fit as the center of the Gaussian and the error bars as the 95% confidence bounds. The q2 feature is fit to two separate Gaussians due to its odd shape (which has no reason to be a simple Gaussian). The q value is extracted as the average of the center of the two Gaussians with the error bars determined by error=√( Gaussian is taken as the 95% confidence bounds. There is no reason to expect such a feature to conform to a single equation fit, such as a single Gaussian, indeed when comparing our data to that seen in reference [3], one would not expect a single Gaussian to be sufficient. Despite using two separate Gaussians to fit the q2 feature, we feel comfortable treating it as a single feature, due to the presence of only two total features along the same direction in reference [3].

Dispersive Measurements of the 1% Mn Sample:
6. Non-dispersive behavior of the Charge Order and AFM peaks: Figure S8 (b),(e) show line cuts taken from the center through the Q * peaks in the FTs of a 5% and 10% sample respectively. All data extracted from the 5% and 10% samples was taken from the raw data with the exception of the 5% -100meV Bias Q* data which was taken from the symmetrized data shown in figure 3. Each line is a different energy from 200mV to -40mV in 20mV intervals for the 5% (The data here were obtained from a conductance map with a setpoint bias of 200meV) and 100mV to -100mV in 10mV intervals for the 10%. The values below -40mV have been excluded as the center inhomogeneity becomes too large and covers up the Q* peaks. The blue circles are the raw line cut data, while the red lines correspond to a fit using a Gaussian function (for the Q* peak) and a 3 rd order polynomial (as background). (Note: a small peak is present in the 10mV, 0mV, and -10mV cuts for the 10% but its intensity compared to the background is too low for it to be seen clearly. For this reason we have chosen to not include these energies in fig 3(h), as a q value cannot be reliably extracted.) Figure S8 (c),(f) show the raw data and fit for the average of all the energies shown in b,e as well as the raw data and fit for the 40mV in the 5% and the 80mV in the 10%. The green bar is the width of the corresponding Bragg peak for scale. (Note: The FT used for the 5% is taken from a much smaller real space map compared to the 10%. Smaller real space maps yield lower resolution FT and explains why all the peaks in the 5% are much broader even when compared to other 5% peaks such as those seen in the FT in Fig. 3. This map was chosen as the Q* peak is measured over a larger range of energies.) The q value is extracted from the fit as the center of the Gaussian, the error bars are 95% confidence bounds of the center. Note that there is a percieved shift in the q value for the 10% Mn sample across the Fermi energy, which is due to a discrete change in the background for negative and positive energies as can be seen in Fig S8(e). ,(e) Line cuts taken through the Q* peak in the 5% and 10% samples respectively. The blue circles are the raw data for each marked energy while the red line shows a 3 rd order polynomial plus a Gaussian fit to the data. Lines are shifted upward for clarity. (c),(f) Raw data and fit for the average of all energies and 40mV data for the 5%, and 80mV data for the 10%. The green bar is the full width of the corresponding Bragg peak taken as the FWHM of a Gaussian fit to the Bragg peaks of the corresponding data. Figure S9 shows additional line cuts taken along the Q* direction showing Q* for the 5% data shown in Fig 3, with their corresponding fits shown as a dashed red line. The data here were obtained from a conductance map with a setpoint bias of -100meV. The fits and extracted values are done using the same method described above for the Q* peak. The values extracted were used for the -100meV values in Fig 3(h). Figure S10 shows similar line cuts taken along the (π/a, π/a) direction showing qAFM for the 5% data shown in Fig 3, with their corresponding fits shown as a dashed red line. The fits and extracted values are done using the same method described above for the Q* peak. There may be a dispersive feature present between 0.4 and 0.45 that could be related to q8. This feature which appears to move to higher q with positive energies, along with an enhanced background makes it difficult to extract any values from fits at positive energies.
FIG. S9: Line cuts taken through the Q* peak in the symmetrized 5% sample shown in Fig 3, with the corresponding fits shown as a dashed red line. The -100meV setpoint Bias is responsible for the suppressed intensity at energies around half of the Bias (-60meV, -80meV, with the -40meV excluded for this reason) due to the setpoint effect. Lines are shifted upward for clarity.
FIG. S10: Line cuts taken through the q AFM peak in the 5% sample with the corresponding fits shown as a dashed red line. The dashed black line marks the q=0.5 location. Lines are shifted upward for clarity.

Theoretical insights into the electronic states of Sr 3 (Ru 1-x Mn x ) 2 O 7 :
FIG. S11: Theoretical insights into the electronic states of Sr 3 (Ru 1-x Mn x ) 2 O 7 (a) tight binding calculation of the Fermi surface of undoped Sr3Ru2O7. (b) Tight binding calculation of the Fermi surface of 5% Mn doped Sr3Ru2O7 (c), Evolution of the susceptibility at Q* and QAFM with decreasing the RuO6 octahedral rotation. Decreasing the octahedral rotation in the calculations corresponds to increase doping. (d), Real part of the susceptibility for the parent compound showing enhanced peaks at Q* and a local maximum at QAFM. Figure S12 shows the real space conductance map and its corresponding raw FT of the 10% doped sample used in Fig. 4 of the main text. The Q* peaks corresponding to the charge order (CO) are marked in the FT.

Visualizing the structure and symmetry of Charge ordering patterns:
FIG. S12 Visualizing the structure and symmetry of Charge ordering patterns Left: Conductance map for the 10% Mn doped sample shown in Fig. 4 of the main text. The long wavelength inhomogeneity in the map is mostly due to the large amount of Mn dopants. Right: Corresponding un-symmetrized FT of the conductance map with the Q* peaks circled.
The FT reveal a broad intensity peak at its center originating from the large amount of inhomogeneity of the real space, mostly due to the large amount of dopants.
To view the CO patterns in real space we multiply the FT in Fig. S12 by broad (with Full width corresponding to the white circles in Fig. S12) Gaussians centered at the four Q* positions. Taking the inverse FT of this reveals the structure of the CO in real space. Figure S13 (a)-(c) shows the results of this method for all the Q* peaks, only the peaks along the Y direction, and only the peaks along the X direction, respectively. These maps reveal a glassy modulation along both directions. The yellow circle in Fig. S13(a) is the expected CO patch size based on the Q* peak in the FT of Fig. S13. The width of the peak is measured as the full width at half maximum (FWHM) of a Gaussian fit to a line cut through the Q* peak this value is than converted into Å (/FWHM) to find the diameter of the circle.