Algebraic correlations and anomalous fluctuations in ordered flocks of Janus particles fueled by an AC electric field

We study the polar collective dynamics of Janus colloidal particles fueled by an AC electric field. When the density is high enough, the polar interactions between the particles induce a polar orientationally ordered state which exhibits features reminiscent of the Vicsek model such as true long-range order and giant number fluctuations. Independent measurements of the polarity and velocity at the single particle level allowed us to investigate the single particle dynamics within the ordered state. We discovered theoretically-unaddressed statistical properties of the ordered state such as the asymmetric relation of polarity and velocity, enhanced rotational diffusion stronger than in the disordered state, and an algebraic auto-correlation of the polarity. Our experimental findings, at the crossroad of the Vicsek physics and the Active Brownian Particles physics, shed light on the so-far-unexplored physics arising from the interplay between the polarity and the velocity.


I. INTRODUCTION
Active matter systems consume and dissipate energy at the level of its local units to generate systematic motion based on their internal degrees of freedom (e.g. polarities) [1][2][3][4]. After decades of extensive research, it is now known that active matter systems can exhibit a wide range of phenomena such as phase separation [5][6][7][8], active chain formation [9,10], active turbulence [11][12][13][14], etc. Among these topics, orientational order in active matter has gained considerable attention since the introduction of the Vicsek model in 1995, where pointwise polar particles move at constant speed toward their polarity and try to align their polarities with local neighbors under the presence of noise [15]. Theoretical work concerning the Vicsek-style models has been done including particle, kinetic and hydrodynamic levels [16][17][18][19][20][21][22][23][24]. In particular, the phenomenological hydrodynamic equations by Toner, Tu, Ramawsamy and their coworkers have predicted the existence of true long-range order and giant number fluctuations (GNF) for the orientationally ordered state [16-18, 23, 25]. In these theoretical and numerical studies, each particle's polarity is assumed to be parallel to its velocity. In experimental systems, on the other hand, particles do not necessarily move toward their polarities due to thermal and/or athermal noise and inter-particle (hydrodynamic, excluded volume, or electrostatic) interactions (e.g. systems of elongated bacteria [26], vibrated polar disks [27], motility assay [28], and active rolling colloids [29,30]).
One of the minimal models for investigating such active systems whose particle velocities do not necessar- ily coincide with their polarities is the Active Brownian Particles (ABP) model. In the ABP model, particles self-propel with a constant self-propulsive force under rotational/translational noise. The particles interact only by volume exclusion and do not exhibit any aligning interactions [31,32]. In spite of the absence of alignment, intriguing collective phenomena such as motility-induced phase separation have been reported for ABP [33]. However, it still remains unclear what lies between non-aligning ABP-like systems and the aligning Vicsek-style models. Active colloidal systems similar to the ABP model have been experimentally realized and studied by using e.g. auto-catalytic colloidal particles in hydrogen peroxide and rolling colloids under a DC electric field. However, active particles fueled by chemicals do not exhibit ordered states reminiscent of those in Vicsek-style models so far [6,7,34]. In the case of rolling colloids, because the particles do not possess intrinsic polarity and that their polarities are only electronic and invisible, in-depth study of the interplay between the polarity and velocity in the ordered state is not possi-ble [30]. Thus, our understanding of how polarity and velocity interact with each other and how they behave in orientationally long-range ordered colloidal systems with hydrodynamic interactions remains elusive. For example, in the systems of elongated bacteria and motility assay, although true long-range order has been observed, particles could only be detected in a coarse grained manner [26,28]. Janus particles [35], which are colloidal particles with two distinct hemispheres, offer a great platform for investigating active matter experimentally [6,8,9,32,[36][37][38][39]. Among these Janus particles, the ones fueled by an AC electric field are especially suitable for analysis from the view point of statistical physics because they constantly convert energy into motion without running out of energy sources and thus, allow long enough measurements for statistical analysis [8,9,32,37,39,40], which is not the case for auto-catalytic particles due to the consumption of chemical fuels. Under an AC electric field, Janus particles display two different mechanisms of swimming depending on the field frequency: induced charge electrophoresis (ICEP) at low frequency [41][42][43]; selfdielectrophoresis (sDEP) at high frequency [44,45]. Importantly, the two swimming mechanisms induce motion from the dielectric hemisphere to the metallic hemisphere [8,9,32,37,38,40,43] or vice versa [9,38,40,45], which allows us to define and observe each particle's polarity based on the location of its metallic hemisphere [8,9]. Not only the swimming mechanisms but also the interactions between particles depend on the field frequency and ion concentration, which lead to various types of collective behavior such as a turbulent phase [37], a swarming state with vortices [38], flagella-like movement [9] and active clustering [8]. Thus, Janus particles provide us an ideal venue for experimentally studying the collective dynamics of ordered systems while keeping track of polarity at the single particle level.
Here in this paper, we perform experiments on the quasi-two dimensional system of active Janus particles fueled by an AC electric field. We specifically use the AC field frequency regime where Janus particles swim by sDEP and exhibit a flocking phase. Through the extraction of both the polarity and velocity for each particle, we show that this flocking phase of Janus particles can exhibit true long-range orientational order. This ordered state of Janus particles also exhibited statistical features such as GNF and algebraic correlations for orientation fluctuations which are reminiscent of the predictions of Toner and Tu [16,17,23]. We further investigate the space-time correlations of density fluctuations leading to an estimate of dynamical and anisotropy exponents for the ordered state of Janus particles. In addition, through single particle measurements, we show that the particles counterintuitively exhibit larger diffusion in the ordered state than that in the disordered state for short time scales. We suggest that this anomalous enhanced diffusion stems from the coupling between particle polarity and velocity through computations of the cross-correlation function, which is a dismissed feature in Vicsek-style models. Our work, through the observation of both velocities and polarities, unravels hitherto unexplored statistical properties within orientationally ordered active colloidal systems.

II. EXPERIMENTAL SYSTEM
Our experimental system is composed of silica based Janus particles with 3.17 µm diameter [9]. To create distinct hemispheres, 35 nm of titanium followed by 15 nm of silica were deposited to a monolayer of silica colloidal particles (Bangs Laboratories, Inc., SS05N) using electron beam deposition and thermal evaporation, respectively. The additional 15 nm layer of silica was deposited to suppress oxidization and the subsequent changes of electrical properties and color (from black to transparent) of the Ti surface. After deposition, the particles were suspended into deionized water and washed by sonication to prevent particles from forming clusters. This washing procedure was repeated three times after waiting for 30 -60 minutes for particle sedimentation. The washed Janus particles were suspended in deionized water and then sandwiched by two 25-nm-silica-deposited indium tin oxide (ITO) coated electrodes (Mitsuru Optical Co. Ltd.), separated by 130 µm thick double sided tape. Before applying the AC electric field, we waited for a few minutes to ensure that the particles have sedimented to the surface of the bottom electrode, which leads to a quasi-two dimensional system of Janus particles ( Fig. 1). Accordingly, the induced-charge electroosmotic flow caused by the electric field lets the particle's polarity turn to the horizontal xy-plane, leading to confined motion in the horizontal 2D-plane [47]. A sinusoidal AC electric field with f = 1 MHz frequency and 16 V pp voltage was applied for all experiments, where V pp denotes the peak-to-peak voltage. In this frequency regime, it is known that Janus particles tend to move toward its metallic hemisphere, which contradicts with the ICEP theory [41][42][43]. This was first observed in [40], and is now understood as sDEP [44,45]. It should be noted that the ordinary electrophoretic effects caused by the negative charge of the silica particles could be negligible in our system since we use an AC electric field. Therefore, ordinary electrophoresis alone cannot compete with the ICEP/sDEP mechanism. Experiments were performed under two different density conditions: the low density regime where ρ 0.78 particles/1000 µm 2 , and the high density regime where ρ 13 particles/1000 µm 2 .
The observations of Janus particles were done using an inverted microscope (Nikon ECLIPSE TE2000-U) with objective lens ((Nikon Plan Fluor ELWD, 20×, NA=0.45 and 40×, NA=0.60 for the low and high density regime, respectively). A halogen lamp was used to illuminate the particles and a green filter was inserted in the light path before the sample to increase the contrast of the obtained images. Under this condition, the polarity of each Janus particle can be distinguished using its dark-gray hemisphere indicating the titanium side, and the light-gray hemisphere indicating the silica-only side (for example see Fig. 3(a), Supplemental Movie 2, 3). The dynamics were recorded at 10 fps and 15 fps for the low and high density regime, respectively, using a CMOS camera (Baumer LXG80, 3000 × 2400 pixels 430µm × 340µm (40×)). This allowed us to capture the large scale statistics of the particles while detecting the polarities and velocities of each particle at the same time.
Image analysis was performed using Trackpy [48] and ImageJ [49]. Specifically, the polarities of the particles were detected by taking the difference of the centroid and center of mass of image intensity for each particle using the built-in algorithms in ImageJ [46]. The duration of the analyzed movies were 360 s and 212 s for the low and high density regime, respectively.

A. Binary polar interaction and following motion
To probe the interaction of the Janus particles in our system, we first measured the relative polarity distribution P (∆θ) for all Janus particle pairs in the low density regime ( Fig. 2(a)). To examine the local interactions of the particles, particle pairs within 6 µm and 30 µm distance were extracted, respectively, and the relative polarity angle ∆θ was measured for each pair ( Fig. 3(a) . The results show that particles within 6 µm range had a narrower distribution with a peak near ∆θ = 0, which suggests the existence of polar interaction within the Janus particles. Indeed, when observing two Janus particles within a short distance, following motion where one particle tends to follow the other particle could often be observed (for example, see Fig. 3(c)). The observed polar interaction and following motion can be understood as a result of the balance of electrostatic and hydrodynamic interactions between the particles. Although the interaction mechanism of Janus particles under an AC electric field is still under discussion [9,38], our previous study experimentally validated the existence of polar interaction induced by the effective electrostatic dipole-dipole interaction between two particles both in the ICEP and the sDEP regimes [9]. Note that, although our study in the sDEP regime uses a similar AC field frequency as in [9], the smaller ion concentration in this experiment than in our previous study generates a longer screening length (Debye length) λ D ∝ 1/ √ c, where c is the ion concentration [9]. The slip velocity u s generated in the electric double layer scales as u s ∝ λ D , following from the Helmholtz-Smoluchowski formula [41]. This results in a weak but relatively stronger pusher-type flow-field than in our previous study [9] which prevents the particles from forming active chains in the current setup. It should also be noted that flow fields calculated based on pusher-type flow fields also show following motion due to the reorientation of the following particle, which is consistent with our observations [50,51]. This synergy of electrostatic and hydrodynamic interactions leads to polar alignment, allowing us to observe the ordered dynamics of Janus particles as we see below.
B. True long-range order Janus particles exhibited distinct global patterns when we changed particle density. While Janus particles exhibited uncorrelated disordered motion in the low density regime ( Fig. 2(a), (c), Supplemental Movie 1), orientationally-ordered coherent motion was observed in the high density regime ( Fig. 2(b), (d), (e), Supplemental Movie 3, 4) due to the local polar interaction. This density dependency of Janus particles is consistent with the previous experimental results of [38]. Note, if we compare the Janus particle system with Vicsek-style models, the density and the noise amplitude should be the basic factors which control the transition from the disordered state to the ordered state. Although the effective noise amplitude, namely the Péclet number, could be adjusted by changing the particle's velocity via the AC electric field voltage, we have not been able to investigate this in the current study. This is because the particles easily adhere to the surface of the ITO electrode at higher voltages forbidding long time observation. Ideally, the effective noise amplitude could also be changed by using particles with different diameters, but it would be nearly impossible to continuously change the noise amplitude in this case. It should also be noted that the global orientation of the system varied each time we conducted the experiment, and thus, is not deterministic. This is analogous to the setting of the Vicsek-style models where the global orientation results from spontaneous symmetry breaking of the continuous rotational symmetry of the system.
To quantitatively assess the properties of orientational order in the low and high density regimes of the system, we computed the polar order parameter, for different sizes of rectangular (5:4) regions of interest (ROIs) taken in spatially-uncorrelated locations (see section SII in the supplemental material for details), where θ is the particle polarity and · S and · t denotes the average within the same ROI and statistically independent time frames, respectively ( Fig. 4, [46]). In the low density regime, the polar order parameter decays as a power law C ∝ 1/ √ S ( Fig. 4(a)). This is consistent with the case with N randomly oriented particles. In contrast, in the high density regime, we can observe a decay slower than a power law within S ∼ 10 5 µm 2 . This asymptotic behavior characterizes the true long-range order, which was previously reported in Vicsek-style models [4] and in the bacterial experiment [26] where the orientational order parameter converges to a positive finite value. Note, the decay slower than a power law in our Janus particle system is observed for about 2.8 decades, which is larger than the range reported in the long-range ordered phase of bacteria [26] and even comparable to the range reported in numerical simulations (see Fig. 4a of [4]). To confirm the existence of true long-range order, we performed nonlinear fitting for the polar order parameter at S < 10 5 µm 2 , using the curve fit function from the scipy package [52]. The data was in good fit with an algebraic convergence to a finite value C −C ∞ ∼ S −γ/2 where C ∞ = 0.74 and γ/2 = 0.53 (see the inset of Fig. 4(b)), which indicates true long-range order in the system. It should be noted that recent simulations showed γ = 0.64 for the Vicsek model [4]. At the same time, however, the decay of the order parameter deviates from this trend of algebraic convergence at S > 10 5 µm. This is caused by the existence of counter-flows of the particles (see supplemental movie 3, 4) which is unavoidable due to the closed boundary conditions and the number conservation of the particles unless we elaborate an effective periodic boundary condition by e.g. a racetrack-like confinement [29]. However, in the case of our electrically driven Janus particles, micro-structures fabricated by e.g. photoresist distort the electric field, which leads to unwanted complex 3D electroconvection and destroys the two-dimensionality of the system. Therefore, we opted for using a simple system by sandwiching as large a droplet of the suspension of Janus particles as possible between the two electrodes. Overall, the finite size scaling analysis assures that the Janus particles exhibit polar collective motion with true long-range order.

C. Interplay of polarity and velocity
To investigate the detailed dynamics of single particle trajectories, we computed the polarity auto-correlation function f p (τ ) for the disordered and ordered states, respectively. Here, where n p (t) denotes the polarity for the particle at time t, and · t takes the average for particle trajectories longer than 24 seconds (Fig. 5(a)(b)). We can first observe that the polarity auto-correlation decays exponentially in the disordered state. At small time scales, this decay was slightly faster than the theoretical estimation for an Active Brownian particle f p (τ ) = exp (−2D r τ ) [31,36], where 1/D r ∼ 24.2 s is the rotational diffusion constant ( Fig. 5(a)) from the Einstein-Stokes relation, which is due to the polarity detection errors [32]. We thus focus on the time scale τ = 0-12 s, where the decay of the auto-correlation exhibited a slope of 21.15 ± 0.08 s, close to that of the theoretical estimation.
In the ordered state, on the other hand, f p (τ ) exhibited a slow power-law-like decay in the long time scale (τ > 1 s, Fig. S4). Since particle polarities in the ordered state tend to align to the global order, we suspected that this power-law behavior of f p (τ ) is related to the fluctuations of the global order. To assess this hypothesis, we compared the connected polarity auto-correlation function in the ordered state with that of the global order (Fig. 6). In Fig. 5(c), we show the connected polarity (velocity) auto-correlations, where ∆n p,v (t) = n p,v (t) − n p,v (t) , and n p,v (t) is the direction of global order defined by polarity or velocity, respectively. A power-law-like decay with a slope close to −0.7 for f p,C (τ ) could be observed for both the polarity and velocity auto-correlations. Interestingly, the connected auto-correlation of the global order f global,C (τ ) [46] also showed a power-law-like decay with a slope of -0.85, close to that of f p,C (τ ) (Fig. 6(b)). Although there is, to the best of our knowledge, no theory predicting the slope for these power law decays in global order and polarity (velocity) auto-correlations, these observations support our hypothesis that the power-law behavior of the polarity auto-correlations is a consequence of the fluctuations in global order. However, it should be noted that in a finite system, the direction of global order diffuses such that the auto-correlation function decays exponentially [24]. Therefore, we should expect that the powerlaw behavior of f global,C (τ ) that we observe is a transient to an exponential decay. Indeed, when we focus on longer time scales (τ > 10 s), we can observe a decay faster than a power-law for f global,C (τ ) (Fig. S5). Although the limited length scale of our experimental system hinders the long time behavior of single particle trajectories, we expect that the connected polarity auto-correlation function f p,C (τ ) will also show a crossover to an exponential decay, reflecting the fluctuations of the global order.
The polarity auto-correlation function for particles in the ordered state also exhibited peculiar behavior at short time scales (τ ∼ 1 s). Fig. 5(b) shows a magnified view of the behavior of f p (τ ), and it can be observed that f p (τ ) of the ordered state decays faster than that of the disordered state at short time scales. Because this time scale is similar to the time scale of the mean free time (distance to the nearest neighbor ∼ 7 µm)/(average velocityv o ) ∼ 0.8 s, we speculated that this enhanced diffusion was caused by the frequent and continuous interactions among closely distanced particles in the ordered state. To confirm this speculation, we computed the polarity-velocity cross-correlation function, where n v (t) = v(t)/ |v(t)| is the unit velocity vector at time t (Fig. 7). The assessment of g pv (τ ) revealed the asymmetric relation of polarity and normalized velocity in the ordered state, where polarity tended to correlate the most with the velocity of 0.2 second in the past (τ ∼ −0.2). In other words, the polarity in the ordered state tended to follow the direction of the velocity of the particle with a finite time delay. It should be noted that, although not statistically significant, g pv (τ ) in the disordered state showed a weak, if any, asymmetry around τ = 0 as well. These results imply that the asymmetry in g pv (τ ) is amplified in higher density situations. We therefore suspect that this asymmetry is caused by electrostatic and/or hydrodynamic torques induced by collisions between individual particles. Overall, our results suggest that the particle polarities tend to follow the particle's velocity, yet local interactions between nearby particles disturb their alignment causing a time delay especially at the high density state. These effects might explain the enhanced diffusion observed in the ordered state.

D. Giant number fluctuations in the ordered state
Giant number fluctuations (GNF), where the variance of the number of particles ∆N 2 grows faster than its mean N (∆N 2 ∝ N 2α with 2α > 1) in long-range ordered phases, are a hallmark of the Toner-Tu phase [16,17,23]. Since its theoretical prediction by Toner, Tu, Ramaswamy and their coworkers [16,17,25], the presence of GNF has been extensively studied in numerical systems [18,21,24,53]. However, in addition to experimental difficulty in realizing collective motion with true long-range order, experimental detection of GNF requires delicate treatment due to many possible pitfalls such as difficulty in disentangling fluctuations arising from e.g. phase coexistence or boundary effects as critically discussed in [26]. Therefore, GNF as the hallmark of the Toner-Tu-Ramaswamy phase had remained undetected until recent experiments on a long-range ordered nematic phase of elongated bacteria [26] and then on the polar ordered phase of colloidal rollers [30]. Here in our Janus particles system, as we have confirmed the true longrange order, it is natural to question whether the ordered state shows anomalous number fluctuations.
To quantify number fluctuations, the mean N and standard deviations ∆N of the number of particles for different sizes of ROIs were calculated for statistically independent frames (Fig. 8). For the disordered state, normal number fluctuations with ∆N ∝ N 0.5 were observed, which is consistent with the central limit theorem. In contrast, anomalous number fluctuations were observed in the ordered state where ∆N ∝ N α with α > 0.5. In fact, the estimated exponent was α = 0.790(5) with 95% confidence level, which is very close to the theoretical prediction α = 0.80 [16,17] and the recent numerical estimate α = 0.84 for the Vicsek model [24]. It is worth noting that the existence of non-trivial GNF has to be discussed under the presence of true long-range orientational order. Indeed, there are several experimental studies that report GNF for systems that are not in the fully ordered state (for details see the supplementary material of [26]). Here, we have set all ROIs to be placed at the center of the whole image so that the ROIs correspond to those that were used for the computation of the order parameter. Accordingly, the gray dashed line in Fig. 8 indicates the length scale S ∼ 10 5 µm 2 at which the order parameter in the ordered phase starts to deviate from algebraic convergence to a finite value (i.e. true long-range order). Therefore, the observation of anomalous fluctuations below this line imply that ordered Janus particles do exhibit GNF in the context of the Toner-Tu phase.

E. Algebraic scars of long-range order
When a continuous symmetry is spontaneously broken, the entire ordered phase should also be characterized by an algebraic decay of its connected correlation functions [54]. In the case of the Vicsek model, a power-law decay for the orientation fluctuation correlation function exemplifies the spontaneous symmetry breaking of the rotational symmetry of the system, and its exponents have been predicted by the hydrodynamic theory of Toner and Tu [16,17,23]. Interestingly, a recent numerical study of the Vicsek model suggested that the ordered state could exhibit GNF exponents close to the conjectured values even though there existed a discrepancy in the exponents of the correlation functions [24]. Thus, in addition to number fluctuations, the fluctuation correlation functions are important measures to probe the statistical properties of orientational order in the system.
To address the questions of whether the orientational fluctuations of Janus particles exhibit scale-free correlations, and if so, whether they follow exponents of the Toner-Tu conjecture, we measured the orientation fluctuation correlation function, where δn p⊥ (t, r) = n p (t, r) − n p is the deviation of a single particle polarity n p from the global order n p , and · r takes an average for all |R| distanced particle pairs (see Fig. 9). Here, the global order n p was calculated by averaging the polarity n p both spatially and temporally. Figure 9(b) shows the decay of C C (R) in the transverse and longitudinal directions to global order. Although there is an anisotropy in the absolute value of the correlation function, it decays algebraically in both directions with the same exponent close to −0.52. This algebraic decay signifies the spontaneous symmetry breaking of the rotational symmetry of the system and consequent scale-free features. However, the exponents of the decay differ from the Toner-Tu conjecture which predicts an anisotropic decay of C C (R) ∝ R 2χ ⊥ , R 2χ/ξ , where χ = −0.2 and ξ = 0.6 are the roughness and anisotropy exponent, respectively [16,17,23]. This discrepancy might stem from the difference in microscopic properties including steric and hydrodynamic interactions as we will discuss in the following section. Interestingly, recent experimental [30] and numerical [24] studies for polar orientational order have both reported weak, possibly vanishing, anisotropy for the orientational fluctuation correlation function. In addition, our analysis on the space-time correlations of density fluctuations also suggests weakly anisotropic sound modes propagating in our Janus colloidal flocks (Fig. S2). This leads to the estimate of the dynamical exponent z ∼ 1.2 while the anisotropy exponent ξ ∼ 1, suggesting a non-diffusive but weakly anisotropic dependence that qualitatively agrees with recent numerical simulations [24]. It should also be noted that the roughness exponent χ also takes similar but distinct exponents for theory, simulation and experiments where χ = −0.20 for the Toner-Tu conjecture [16,17,23], χ ∼ −0.31 for numerical simulations [24], χ ∼ −0.38 for active Quincke rollers [30], and χ ∼ −0.26 for our system of Janus particles. These differences in the exponents might indicate that some of the nonlinearities that were neglected in the original calculations [16,17] are relevant asymptotically, as suggested in [24]. Our measurements of C C (R) and the sound modes, along with recent numerical [24] and experimental [30] studies, may pave the way for further theoretical work on collective systems.

IV. DISCUSSION
We studied the collective motion of active Janus particles fueled by an AC electric field. At high density, the Janus particles exhibited polar orientational order which allowed us to investigate and compare the statistical properties of the system with that of the Toner-Tu phase. We have confirmed GNF within the polar orientational ordered state with an exponent similar to but different from the theoretical predictions of Toner and Tu. In addition, an algebraic decay was observed for the orientation fluctuation correlation function which is also associated with the spontaneous symmetry breaking of the rotational symmetry of the system. The observed exponents showed a discrepancy with the theoretical predictions and recent numerical simulations [24] of the Vicsek model, suggesting a different mechanism for the emergence of long-range order. As pointed out in [26] and then both experimentally and numerically verified in [28,55], weak volume exclusion that allows crossing over other particles is considered a key to make collective motion of elongated objects similar to the Vicsek-style interactions and thus achieve long-range order at least for the otherwise turbulent or clustering 2D nematic systems such as bacteria and gliding microtubules [26,28]. In the case of polar collective motion, the weak exclusion may not be a requirement for the long-range order because the particles cannot overlap each other in our experimental setup as in granular experiments [27,56]. What realizes long-range order in polar systems could be the following motion exemplified in Fig. 3, as discussed previously in e.g. mutant Dictyostelium cells [57], and/or the torque aligning the polarity to the velocity, as studied in vibrated polar disks [58].
Our results might suggest that, although GNF could be observed, algebraic decay as in the Toner-Tu phases are difficult to observe in experimental systems as well as in numerical simulations. The exponents of the algebraic decay for the orientation fluctuation correlation function differed from that observed for polar flocks of Quincke colloidal rollers [30]. Although Janus particles and Quincke rollers are both considered as polar particles with polar interaction, the polarity acquisition mechanisms are completely different. Whereas the polarity of Janus particles is a priori defined by its metallic hemisphere, Quincke rollers acquire polarity by the spontaneous symmetry breaking of charge distribution on the particle scale. It could be interesting to investigate whether this difference in the polarity acquisition mechanism alters the global properties of the orientationally ordered system because the latter might have additional soft modes on the particle level.
Tracking both the velocity and polarity of all particles enabled the detailed assessment of the single particle dynamics both in the disordered and ordered states. Strikingly, the polarity auto-correlation showed a slow algebraic decay in the polar ordered state, which has never been studied even numerically. It could be interesting if this slow decay in polarity could be observed not only in active colloidal systems but also in collective biological matter [57].
Our results also suggest the interplay between the polarity and velocity of the particles which might be explained by the hydrodynamic and/or the electrostatic torques between the particles. We have shown that this interplay is amplified in the high density ordered state. This interplay between polarity and velocity is something that cannot be observed in the polar Vicsek model where the two vectors are identical and degenerated. A recent numerical study shows that the global properties could differ drastically by decoupling the polarity and veloc-ity of the system [59]. Our polarity resolved realization of the collective phase of self-propelled Janus particles can contribute to bridging the understandings on the two major active matter classes: the Vicsek physics and the Active Brownian Physics. Overall, our system of Janus particles provides a large scale statistical study for polar ordered active colloidal particles, which opens a new avenue on further theoretical and experimental works for the field of active matter. To investigate the sound modes in the ordered state of Janus particles, we computed the power spectra using the particle coordinates [30]. The spatial Fourier transform of the density field is defined as, where x j (t), y j (t) are the instantaneous particle coordinates, j takes an average of all particles and q(cos θ, sin θ) is the wave vector making an angle θ with the global order θ . Here, the global order is defined as the mean polarity of all particles averaged through all frames. The power spectra |ρ q,ω | 2 (Fig. S1) were obtained by performing time Fourier transformations to the two-time autocorrelations of the density fields ρ q (t) · ρ † q (t + τ ) t . The speed of sound modes c(θ) = lim q→0 [ω (θ) /q] were obtained from the positions of the peaks in the dispersion relations given by the power spectra (Fig. S2). Here, the positions and widths ∆ω of the peaks as a function of q in the power spectra were determined by fitting the |ρ q (ω)| 2 curve using a Lorenztian function, The θ dependency of c(θ) could be fitted by a cosine curve: a cos (bθ + δ) with a = 9.5, b = 1.0, and δ = 0.65, respectively ( Fig. S2(a)). The speed of sound c(θ) is comparable to the mean velocity of the particles in the ordered state (v o ∼ 8.5µm/s) which indicates that the observed density fluctuations are caused by the advection of the particles. In addition, sound modes of particles that move in the same direction with a uniform velocity exhibit a sinusoidal dependency against θ, which again suggests that the sound modes observed in our system originate from the advection in the system.
We have also plotted the peak widths ∆ω, characterized as the full width at half maximum (FWHM) of the fitted Lorenztian, as a function of q (Fig. S2(b), (c)). Although the statistical noise in our system kept us from observing the two mixed sound modes as predicted by Toner and Tu, we were still able to observe the algebraic behaviors underlying the sound modes. For both longitudinal and transverse directions, ∆ω(q) showed a slope close to 2 where 3 µm < 1/|q| < 7 µm, indicating that the sound modes are subject to diffusive damping. On the other hand, ∆ω(q) showed a slope smaller than 2 at length scales larger than typical interparticle lengths 7 µm < 1/|q| which suggests the existence of long-range correlations that suppress diffusive damping. Interestingly, ∆ω(q) showed similar exponents for both longitudinal and transverse directions which is consistent with the results of C C (R) (Fig. 9), indicating weak anisotropic correlations with respect to global order. Note, the theoretical predictions given by Toner and Tu are ∆ω ∼ q z/ξ , q z ⊥ (q ,⊥ → 0), where z = 1.2, ξ = 0.6 are the dynamical and anisotropy exponent, respectively. The red, blue and yellow areas correspond to length scales smaller than the particle's diameter (< 3 µm), length scales larger than the typical interparticle length in the ordered state (> 7 µm), and length scales in between, respectively.

SII. CALCULATION OF THE POLAR ORDER PARAMETER
To assess the relation between C and S when the center of ROI is taken in different places compared to that when the ROI is fixed to the center of the whole field, we have calculated C with ROIs fixed to the center, top left, bottom left, bottom right, and the top right of the whole field, respectively ( Fig. S3(a)(b)). This analysis indicated that using ROIs only from a certain location gives systematically erroneous results for C. Therefore for every ROI size, we have divided the whole view into grids and sampled ROIs from various locations to reduce systematical errors for C (Fig. S3(c)). Here, we have sampled a ROI in every three grids (both in the horizontal and longitudinal direction) to prevent spatial correlations among the sampled ROIs. This results in a total number of 36 ROIs to average, for example, when S ∼ 500 µm 2 . The frames were also sampled every 100 frames (∼ every six seconds in the ordered state) to prevent temporal correlations. The results for C are shown in Fig. 4. C based on ROIs sampled from the whole view shows an algebraic decay to an asymptotic value C ∞ = 0.74 with a slope of S −0.53 (Fig. 4(b)). The result where C based on ROIs sampled from the whole view exhibits a slower decay than a power law, indicates that the asymptotic behavior of C is robust.  Fig. 4. We have divided the whole view to grids and sampled a ROI in every three grids (both in the horizontal and longitudinal direction) to prevent spatial correlations within the sampled ROIs.

SIII. AUTO-CORRELATION FUNCTIONS
The polarity auto-correlation function f p (τ ) in a log-log plot is shown in Fig. S4. A slow power-law-like decay for τ > 1 s with a slope close to −0.15 can be observed for f p (τ ) in the ordered state. The auto-correlation function of the polarity fluctuations from global order f global,C (τ ) = ∆ñ p (t)·∆ñ p (t+τ ) t , where ∆ñ p (t) = n p − n p t,ordered ordered , is shown in Fig. S5. f global,C (τ ) shows an oscillation with a period of ∼ 20 s while it decays. This is consistent with the time scale of the time needed for the particles to propel through the ordered areav o × 20 s ∼ 170 µm which is close to the length scale where the order parameter deviates from true long-range order.

SIV. POLARITY DETECTION VIA IMAGEJ
In order to obtain the polarity of each particle, we have applied image analysis methods using ImageJ [49]. Here we give a description of the whole image analysis procedure in order to obtain the orientation which is defined by the metal hemisphere. The main idea is to binarize the images of Janus particles and obtain the orientation of each particle by taking the difference of the center of mass (intensity) and the centroid. In Fig. S6, we show how the processed image will look like when we perform the following procedure on the raw experimental image (Fig. S6(a)).
(a) (b) (c) (d) (e) (f) (g) (h) frames for computing the velocity distribution in Fig. 2(c)(d). For the computation of the auto-correlation function f v (τ ), the velocity was exceptionally defined as, v(t) = x ci (t+3)−x ci (t), in order to reduce noise in the measurements by using a larger window. It should be noted that, we used the symmetric definition v(t) = x ci (t + 1) − x ci (t − 1) for the computation of the cross-correlation function g p,v (τ ) to avoid confusion in the interpretation. The order parameter in Fig. 4 was computed using the particle polarities detected by ImageJ. The time series of the total number of particles in the whole field of view and the polar order parameter for the whole view are given in Fig. S7.
For inferring binary interactions in Fig. 3, the locations of particle pairs within 6 µm and 30 µm distance were extracted for statistically independent frames which were 25 frames (= 2.5 s) apart. Janus particles in the orientationally disordered state. The particles were exposed to a 1000 kHz, 16 V pp AC electric field. The movie shows the first 1800 frames out of the 3600 frames used for binary collision and other statistical analysis. The resolution is scaled down to 640 × 512 from the original 3000 × 2400 to reduce the file size. The movie is played at three times the real speed.
Supplemental Movie 2: (2 Disordered crop.mp4) The same experiment as Supplemental Movie 1 with a cropped view of 562 × 382 pixel view for the 15th -55th frames. The original resolution is kept to ensure the visibility of each particle's polarity. The movie is played at the real speed.
Supplemental Movie 3: (3 Ordered.mp4) Janus particles in the orientationally ordered state. The particles were exposed to a 1000 kHz, 16 V pp AC electric field. The movie shows the first 1800 frames out of the 3180 frames used for statistical analysis. The resolution is scaled down to 640 × 512 from the original 3000 × 2400 to reduce the file size. The movie is played at two times the real speed.
Supplemental Movie 4: (4 Ordered colored.mp4) Janus particles in the orientationally ordered state with colored by their deviation from global order. The movie shows the 1500th -3000th frames out of the 3180 frames used for statistical analysis. The colors show the polarity deviation for each particle from the global order in radian. The movie is played at four times the real speed.