Dynamic evolution of hyperuniformity in a driven dissipative colloidal system

The experimental system is the same as in references [41, 42], where detailed explanations are provided on the setup and its operational dynamics. As can be seen from figure 1(a), the setup consists of a femtosecond laser (Spectra-Physics, Spirit One 1040-8-SHG, 1040 nm wavelength, 1 MHz pulse repetition frequency, and ∼6 μm spot size) coupled to an inverted microscope (Nikon, Eclipse Ti-U). The experiments are performed on a quasi-2D-confined, monodispersed pure polystyrene spheres (1 μm in diameter) suspended in water (2 wt%, microparticles GmbH). Experiments are recorded using a CMOS camera (DCC1645C-HQ, Thorlabs) with 1280 × 1024 pixels and a rate of 15 frames per second. A blue LED is used for imaging illumination.

Zoom In Zoom Out Reset image size

As detailed in references [41, 42], initially, the system has a statistically isotropic distribution of Brownian particles. Upon turning on the laser, a cavitation bubble forms along with laser-induced Marangoni flows, which drag the particles towards their aggregation at the bubble boundary (figure 1(b), left). Turning off the laser stops the drag flow, and particles disassemble back into the system, eventually relaxing to an isotropic distribution (figure 1(b), right).

Thirty experiments have been performed for this study. Although we confirmed our central claims for all experimental video recordings, combining these experiments to provide ensemble averaging of specific quantities such as number variance was not possible because each experiment is necessarily different as neither the sequence nor the exact location of lattice imperfections or external perturbations can be repeated in this strongly stochastic system. Nonetheless, the evolution of hyperuniformity from class I to class II, class III, and non-hyperuniform states was observed and quantified for each experiment. Due to the sheer amount of data, we cannot provide the results of all experiments here. Therefore, we present our analyses over one experimental video, which can be viewed from video 1 (https://stacks.iop.org/JPCM/33/304002/mmedia). We further offer appendix figure 11 and appendix table 3 that show the same analysis as in figure 3 performed on a different experiment.

Video 1 is divided into four parts since camera recording stops automatically upon hitting the file size limit of 4 GB. The original recording of ∼1 h is sped up 15 times to reduce the duration to roughly 4 min. The experimentalist immediately starts recording a new video with roughly 20 s of data acquisition gaps between each part.

The first part of the video shows the formation of a nearly perfect hexagonal crystal at a cavitation bubble boundary. The red square frame in the video (orange square frame in figure 2(b)) indicates the field of view (720 × 720 pixels) where we perform our analyses. The analyses start at t = 6.5 min (frame # 5850) (figure 1(b), left), and end at t = 68.3 min (frame # 61425) (figure 1(b), right). The laser position is denoted in figure 1(b), left.

Zoom In Zoom Out Reset image size

As can be seen from the video, the system is highly dynamic. Particles never stop their motion even when they are embedded into the crystal. In part 2, the experimentalist turns off the laser, and the crystal disassembles. As will be discussed in section 3, during this disassembly, various lattice imperfections such as line and point defects, dislocations, and voids emerge, propagate, and annihilate because physical boundaries change (e.g., bubbles shrink and vanish), additional vectorial fluid flows emerge, and the particles' motions increase. All these changes are reflected in the spatiotemporal distribution of the particles and the interparticle correlations.

Upon turning off the laser, crystal disassembly starts from the lower corners of the field of view (figure 4(b)) resulting an inhomogeneous spatial and velocity distribution of particles. Particles at the lower corners move faster than those in the middle and at the upper corners. Roughly at around t = 18 min (frame # 16229), half of the field of view is solid-like and the other half is liquid-like (figure 4(c)). Our analyses continue in parts 3 and 4 of video 1 until t = 68.3 min (frame # 61425).

The system has relatively fast dynamics, and to assess their effect on hyperuniformity we developed a computational toolbox to perform real-time analysis on video 1. This toolbox provides the simultaneous characterization of hyperuniformity in real- and reciprocal-spaces together with the evolution of several different order metrics and measurements showing the effect of external perturbations on particle movements. On-the-fly analyses of the experiments can be viewed from video 2, a logic map to these analyses is given in figure 2(a), and appendix figure 12 shows a representative screen view of the toolbox.

2.2.1. Particle detection

2.2.2. Dynamic Lindemann parameter map

The second window at the top row of video 2 and appendix figure 12 is the dynamic Lindemann parameter ${L}_{i}\left(t\right)$ [41, 45, 46], which gives a local measure for interparticle distances through

$$\begin{equation}{L}_{i}\left(t\right)=\frac{\sqrt{{\langle {r}_{i}\left(t\right)}^{2}\rangle -{\langle {r}_{i}\left(t\right)\rangle }^{2}}}{\langle {r}_{i}\left(t\right)\rangle },\end{equation} \tag{ 1 }$$

where $\langle {r}_{i}\left(t\right)\rangle$ is the average distance between particles, i, and its neighbouring particles. The averages can be calculated by $\langle {r}_{i}\left(t\right)\rangle =\frac{{\sum }_{j}^{{N}_{i}\left(t\right)}{r}_{ij}}{{N}_{i}\left(t\right)-1}$, where ${N}_{i}\left(t\right)$ is the total number of particles in the scan area of the particle i (i.e., the circular region with a radius of 2.8 times the particle radius surrounding the particle i), including the particle i itself. j is an integer index scanning over all neighbouring particles.

Values of the dynamic Lindemann parameter per particle lie between 0 and 1, which are also colour-coded. Particles that are distant and independent of each other, similar to the gaseous state, are coloured red, ${L}_{i}\left(t\right)=1$. Particles are coloured blue if their neighbouring particles are at their close-packing arrangements, ${L}_{i}\left(t\right)\sim 0$, as in the solid-state. Particles coloured in between are neither distant and independent of each other, nor they are close-packed, $0{< }{L}_{i}\left(t\right){< }1$, similar to the liquid state or, e.g., particles at the grain boundaries [41]. The colour map of the dynamic Lindemann parameter gives a visualization of the particles' spatial correlations and clearly shows the position and propagation of various lattice imperfections frame-by-frame.

2.2.3. Packing fraction

The last window at the top row of video 2 and appendix figure 12 is the packing fraction, ϕ, which is calculated through Voronoi cell analysis using [47, 48]

$$\begin{equation}\phi =\frac{{N}_{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}}\pi {\left(\frac{D}{2}\right)}^{2}}{A},\end{equation} \tag{ 2 }$$

where N_total is the total number of particles, D is the average particle diameter, and A is the aggregate area calculated using the Voronoi diagram. The calculation of A is performed using Matlab function 'voronoin', then adding up the areas of the inner Voronoi polygons. A plot showing time evolution of the packing fraction can be seen in appendix figure 14.

2.2.4. Particle speed and velocity maps

2.2.5. Fourier and reciprocal space images

2.2.6. Hyperuniformity characterization

The last row of video 2 and appendix figure 12 shows hyperuniformity characterization through local and scaled number variances and static structure factor calculations.

We used monodispersed hard-sphere particles with identical diameters (1 μm on average). Naturally, the system can be treated as a point configuration, and hyperuniformity can be assessed from the number variance [2, 6]. However, we additionally calculate the volume fraction variance since our system comprises colloidal particles suspended in water and can be regarded as a two-phase system. The volume fraction variance analyses are provided in the appendix figure 15, which do not contradict our claims on hyperuniformity based on the number variance analyses.

We calculate the number variance from both real-space and reciprocal-space (i.e., the power spectrum). The two methods are performed simultaneously over the same field of view with a fixed size (720 × 720 pixels) and position throughout the experimental video.

For real-space analysis, we choose half of the field of view (360 × 360 pixels) as the sampling region (figure 2(b)), which is sampled by a series of circular observation windows with varying radii (0.1 ⩽ R ⩽ 10) [6]. We particularly choose circular windows as the hyperuniformity measurements are shown to be affected by the geometry of observation windows and that the circular windows are preferred to simplify the calculations [2, 49].

Centres of the observation windows are randomly distributed within this sampling region (indicated by the red dots in figure 2(b)). The number of windows is fixed to be 2500 for each R value. We use the same centres for all radii for one frame. We also tried number variance calculations with a different number of windows for each R value, the lowest being 50 and the highest 10 000, which did not significantly change our results (data not shown).

Next, we count the number of particles, N, within each observation window as a function of its radius, R, in two-dimensional space. Lastly, for each R value, we calculate the variance of the particle count over all the observation windows, which generates the local number variance ${\sigma }_{N}^{2}\left(R\right)$ as a function of R, where R is in the unit of average particle diameter (1 μm):

$$\begin{equation}{\sigma }_{N}^{2}\left(R\right)\equiv \langle {N}^{2}\left(R\right)\rangle -\langle N{\left(R\right)}^{2}\rangle .\end{equation} \tag{ 3 }$$

For reciprocal-space analysis, we compute the local number variance from the structure factor using [2]

$$\begin{equation}{\sigma }_{N}^{2}\left(R\right)=2\rho \pi {R}^{2}{\int }_{0}^{\infty }S\left(k\right)\frac{{{J}_{1}}^{2}\left(kR\right)}{k}\enspace \mathrm{d}k,\end{equation} \tag{ 4 }$$

where ρ is the number density, k is the wavenumber, $S\left(k\right)$ is the static structure factor, and J₁ is the first-order Bessel function of the first kind. To calculate the scaled number variance, we scaled the local number variance with the window surface area using ${\sigma }_{N}^{2}\left(R\right)/{R}^{2}$.

S(k) is calculated by taking the angular average of the reciprocal-space image. The first two values in the S(k) are set equal to the third value to suppress the forward scattering. Suppression of the forward scattering ($S\left(k\right)$ value tightly close to k = 0) is a well-known procedure [6, 10, 50–52]. The number variance measurement based on $S\left(k\right)$ is already effectively omitting the forward scattering because it is multiplied by the wavenumber vector (4). However, $S\left(k\right)$ is calculated by angularly averaging the 2D reciprocal space that is represented by an image, i.e., a discrete grid with side length L, which is not always odd. Therefore, the forward scattering can show up in the first two values of the $S\left(k\right)$ (at k = 0 and $\left.k=2\pi /\right.L$).

The forward scattering is expected to be expanding for more than 1 pixel when the field of view is limited. The spectral resolution is 1 pixel, so the forward scattering in 1D is expected to influence the first two values within the spectral resolution's error range. We choose to omit forward scattering by smoothing out that very limited central part, thus setting the first and second values of the $S\left(k\right)$ equal to the third value. However, one should note that setting these values equal to zero is wrong as it forces the system to approach hyperuniformity.

In appendix figure 16 we provide numerical simulations in 1D of a periodic lattice, a Fibonacci quasicrystal, the eigenvalues of a single random matrix from a Gaussian unitary ensemble (GUE), and uniform random distribution of points along a line to ensure that such omission does not affect the reciprocal-space number variance measurements. In appendix figure 17, we provide additional numerical simulations for a 2D square lattice. We add uncorrelated perturbations to each point, forming the lattice with varying degrees, along with a uniform random distribution. These numerical simulations show that the omission technique does not affect the agreement of real- and reciprocal-space local number variance measurements in different 1D and 2D cases. These simulations exhibit the correct behaviour for different hyperuniformity states.

Two points have to be clarified for the local number variance graphs shown in video 2 and appendix figure 12: first, along with the experimental number variance curves (black dots for real-space and purple squares for reciprocal-space analyses), we show guide-curves that follow the expected number variance for class I, ${\sigma }_{N}^{2}\left(R\right)\sim R$, class II, ${\sigma }_{N}^{2}\left(R\right)\sim R\enspace \mathrm{ln}\left(R\right)$, and non-hyperuniform, ${\sigma }_{N}^{2}\left(R\right)\sim {R}^{2}$, states.

Second, we are not normalizing our number variance curves shown in video 2 and appendix figure 12, because hyperuniformity is analysed while the crystal is disassembling. Normalization is structure-dependent and requires an order unity. This is not always applicable here since the system switches from an ordered to a disordered state during the experiment. A solution can be found by applying a fitting function to the number variance curves calculated through real- and reciprocal-spaces. Still then, it makes on-the-fly calculations prohibitively difficult and disproportionately increases the processing time for a total of 61425 frames. Therefore, we separately provide the fitting functions for each frame of interest in the graph legends of figures 3, 4, 7, 8 and 10, and in table 1 with R-squared and root-mean squared-error (rmse) values. We also prepared table 2, providing the same information for each frame discussed in the paper but not shown in the figures.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 1. Best-fit function and fitting parameters of the real- and reciprocal-space number variance curves are provided for each frame presented in figures 3–10.

	Time	Frame #	# Particles	Real-space			Reciprocal-space
				Best fitting function	R-squared value	RMSE	Best fitting function	R-squared value	RMSE
Figure 3	6.9 min	6177	2105	0.526R	0.629	1.02	0.570R ln(R)	0.987	0.407
	9.6 min	8607	2027	0.602R	0.844	0.610	0.760R ln(R)	0.984	0.663
	26.5 min	23858	1262	0.523R2	0.991	1.53	0.896R ln(R)	0.991	0.585
	27.7 min	24886	1028	1.06R ln(R)	0.997	0.352	1.31R ln(R)	0.994	0.707
	28.5 min	25617	1016	1.08R	0.975	0.406	1.02R ln(R)	0.996	0.428
	32 min	28767	963	0.755R ln(R)	0.996	0.310	0.988R ln(R)	0.998	0.278
	34 min	30582	910	0.607R1.79	0.994	0.788	1.15R ln(R)	0.998	0.346
	54.8 min	49346	657	0.718R2	0.968	4.16	0.441R1.86	0.999	0.146
Figure 4	9.4 min	8468	2028	0.614R	0.844	0.704	0.665R ln(R)	0.992	0.374
	12 min	10839	1988	0.632R	0.818	0.790	0.392R1.93	0.998	0.386
	16.6 min	14960	1710	0.0281R3.87	0.999	0.871	0.921R2	0.993	0.234
	27.2 min	24478	1045	0.679R ln(R)	0.967	0.750	0.909R ln(R)	0.994	0.439
Figure 5		(a)	—	0.502R	0.471	1.25	0.502R	0.583	1.06
		(b)	—	0.487R	0.471	1.22	0.861R2.57	0.999	3.49
		(c)	—	0.0304R4.72	0.999	5.18	5.46R2.90	0.999	0.382
		(d)	—	0.488R	0.479	1.20	0.585R2.58	0.999	1.08
		(e)	—	No fit	—	—	2.43R2.73	0.999	8.52
		(f)	—	0.380R3.92	0.999	8.96	5.25R3.15	0.999	60.1
Figure 6	(a)	Leftmost	—	0.360R3.27	0.999	4.21	0.620R	0.659	0.814
		Middle	—	0.008 70R4.22	0.999	0.772	1.06R2	0.991	3.23
		Rightmost	—	0.305R	0.555	0.646	0.297R	0.889	0.259
	(b)	Leftmost	—	0.333R3.23	0.999	3.96	0.633R	0.719	0.747
		Middle	—	0.003 95R4.37	0.999	0.611	0.846R2	0.990	2.73
		Rightmost	—	0.342R	0.943	0.206	0.358R	0.962	0.188
	(c)	Leftmost	—	2.01R3.30	0.999	26.0	No fit	—	—
		Middle	—	0.0406R4.32	0.999	2.24	3.14R2.24	0.994	13.8
		Rightmost	—	0.562R2	0.984	2.21	0.353R2	0.998	0.518
Figure 7	10.3 min	9258	2018	0.612R	0.839	0.686	0.292R1.93	0.998	0.327
	10.3 min	9259	2014	0.591R	0.864	0.602	0.827R ln(R)	0.988	0.720
	10.3 min	9265	2015	0.607R	0.838	0.679	0.771R ln(R)	0.986	0.641
Figure 8	8.9 min	7967	2112	0.519R	0.580	1.08	0.540R ln(R)	0.984	0.430
	9.1 min	8006	1936	0.0157R3.5	0.996	0.921	0.426R2.03	0.998	0.630
	9.1 min	8145	2025	0.340R ln(R)	0.894	0.711	0.799R ln(R)	0.988	0.593
Figure 10	28 min	25184	1024	0.903R	0.961	0.492	0.816R ln(R)	0.999	0.169
		(b)	—	0.258R + 0.0626R ln(R)	0.993	0.0842	0.289R + 0.043R ln(R)+0.009R2	0.999	0.0210
		(c)	—	0.840R−0.722	0.988	0.211	0.962R ln(R)	0.989	0.676
		(d)	—	0.824R + 0.118R ln(R)	0.956	0.593	0.696R ln(R)	0.998	0.196
		(e)	—	1.06R ln(R)	0.997	0.359	0.993R ln(R)	0.997	0.380

Table 2. Best-fit function and fitting parameters of the real- and reciprocal-space number variance curves are provided for each frame described in the text, but not shown in the figures.

Time	Frame #	# Particles	Real-space			Reciprocal-space
			Best fitting function	R-squared value	RMSE	Best fitting function	R-squared value	RMSE
19.0 min	17060	1637	0.0355R3.75	0.999	0.790	0.537R2.18	0.997	1.28
20.0 min	18001	1552	0.0360R3.66	0.999	0.966	0.472R2	0.998	0.601
20.7 min	18660	1499	0.0500R3.50	0.999	0.601	0.594R1.8	0.996	0.704
21.3 min	19161	1467	0.0422R3.43	0.999	0.609	0.969R ln(R)	0.989	0.712
25.0 min	22506	1296	0.126R2.58	0.999	0.342	0.585R1.6	0.992	0.625
25.9 min	23350	1275	0.179R2.43	0.999	0.415	0.797R ln(R)	0.998	0.242
35.9 min	32330	898	0.498R1.71	0.991	0.703	0.526R1.7	0.997	0.402
37.6 min	33827	873	0.405R + 0.254R ln(R)	0.997	0.150	0.414R1.8	0.999	0.170
38.0 min	34228	865	0.699R ln(R)	0.988	0.493	1.06R ln(R)	0.998	0.263
43.4 min	39044	788	0.617R1.7	0.997	0.441	0.593R1.75	0.997	0.502
43.7 min	39368	787	0.592R2	0.986	2.12	0.624R1.6	0.996	0.447
43.8 min	39440	786	0.250R2.38	0.999	0.444	0.396R1.94	0.999	0.291
62.6 min	56329	616	0.0206R3.38	0.984	1.69	0.512R2	0.998	0.571
62.8 min	56556	589	0.136R2.87	0.999	0.944	0.579R1.8	0.998	0.412
63.2 min	56848	587	0.139R3.04	0.999	0.806	0.905R ln(R)	0.991	0.587
64.5 min	58056	575	0.0992R3.02	0.999	0.415	0.721R1.87	0.998	0.719
65.7 min	59110	567	0.0618R3.08	0.999	0.708	1.30R ln(R)	0.997	0.499
66.4 min	59725	554	0.0830R2.92	0.994	1.45	0.567R1.88	0.997	0.685
66.9 min	60207	545	0.0519R3.12	0.996	1.23	0.649R1.8	0.998	0.546
68.3 min	61425	526	0.892R1.23(2 ⩽ R ⩽ 6)	0.988	0.200	0.702R1.88	0.997	0.853
—	—	—	0.006 29R3.89(6 ⩽ R ⩽ 10)	0.997	0.625	—	—	—

video 1(frame # 39440)
FOV	# Particles	Length (pixel)		Real-space			Reciprocal-space
			Best fitting function	R-squared value	RMSE	Best fitting function	R-squared value	RMSE
F0	200	360	—	—	—	0.995R	0.954	0.598
F1	207	360	—	—	—	1.25R ln(R)	0.999	0.255
F2	228	360	—	—	—	0.867R ln(R)	0.999	0.0925
	720	680	0.652R2	0.993	1.60	0.432R1.82	0.999	0.146
	786	720	0.250R2.38	0.999	0.444	0.396R1.94	0.999	0.291
	958	800	0.222R2.36	0.999	0.299	0.382R1.92	0.999	0.124
	1360	970	0.416R2	0.994	0.958	0.409R2	0.999	0.402

The fitting functions are generated via the native least-square method-based curve fitting toolbox of Matlab. We assigned a hyperuniformity class to each structure by comparing the R-squared values of different fitting options. When the error values are close to each other (difference < 0.02), we choose the single parameter fit function. Small values of window radii are omitted when fitting the local number variance data, which is applied only to 2 ⩽ R ⩽ 10.

On a final note, while our real- and reciprocal-space analyses indicate disordered hyperuniformity when the crystal is disassembled, the corresponding static structure factor calculations at infinite wavelengths remain nonzero, $S\left(k\to 0\right){ >}0$. In section 3.6 we discuss this issue and conclude that in a finite-sized experimental system with relatively fast dynamics, $S\left(k\right)$, by itself, cannot reliably diagnose the hyperuniformity.

The system starts from a strong hyperuniform state (class I) for an almost perfect hexagonal crystal with only a few lattice imperfections (figure 3(a) and video 2, t = 6.9 min frame # 6177). Corresponding number variance curve shown in figure 3(b) grows linearly as the window radius, ${\sigma }_{N}^{2}\left(R\right)\sim R$, and scaled number variance (figure 3(c)) approaches zero while oscillating. On the other hand, fitting parameters given in table 1 indicate that the number variance curve calculated from the reciprocal-space behaves logarithmically, ${\sigma }_{N}^{2}\left(R\right)\sim R\enspace \mathrm{ln}\left(R\right)$, pointing to a class II hyperuniform state. The reason is discussed in sections 3.2–3.5 with the aid of figures 4–9. We conclude that any lattice imperfections or inhomogeneities that break the system's periodicity will show up in the reciprocal-space calculations. Corresponding Fourier- and reciprocal-space images and the static structure factor confirm the periodic lattice configuration by the presence of intense Bragg peaks (video 2 and appendix figure 18).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Video 2 shows persistent hyperuniformity of the crystal until t = 8.9 min frame # 7967 at which the laser is turned off. Right at this moment, the lattice suddenly relaxes, interparticle distances grow, and several lattice imperfections emerge and propagate, as clearly seen from the dynamic Lindemann parameter map of figure 4(a) and video 2. This is also reflected in the packing fraction graph of video 2 and appendix figure 14 with a sudden drop from ϕ ∼ 0.84 to ∼0.78, which briefly recovers to ∼0.81, but continues to drop since the crystal is disassembling. There are three gaps in this graph, corresponding to the gaps in data acquisition, as discussed in section 2.1.

The underlying cause of the emergent lattice imperfections is strikingly visual in the particle speed and velocity maps. The formation of vectorial flows in various directions results in increased activity in particle displacement and speed profiles. As will be discussed in section 3.5 with the aid of figures 7–9, despite such strong perturbations, the lattice's overall periodicity is preserved, which is confirmed via Fourier image and power spectrum analyses shown in video 2 and appendix figure 18.

As a result of these vectorial flows, hyperuniformity is temporarily lost (video 2 and figures 7 and 8) as the density fluctuations scale with the window surface area and the tail of the scaled number variance curve flattens (converges to a nonzero value) for large window radii (R > 6). The hyperuniformity is recovered shortly (9 < t < 11 min), where real- and reciprocal-space analyses point to class I and class II hyperuniform states, respectively. Representative curves for t = 9.6 min (frame # 8607) are given in figures 3(b) and (c) with their fitting parameters in table 1. Corresponding Fourier- and reciprocal-space images and the structure factor in appendix figure 18.

At around t = 11 min of video 2, it is seen that the crystal disassembles starting first from the lower-left then the lower-right corners, which is also shown in figure 4(b). This is also evident in the particle speed and velocity maps.

The peculiarities in density fluctuations start around t = 12 min. Real-space calculations indicate a class I hyperuniformity, ${\sigma }_{N}^{2}\left(R\right)\sim R$, until t = 14 min whereas reciprocal-space calculations point to a non-hyperuniform state, ${\sigma }_{N}^{2}\left(R\right)\sim {R}^{2}$. At around t = 17 min half of the field of view is solid-like, the other half is liquid-like (video 2 and figure 4(c)). As expected, such a significant inhomogeneity in spatial particle distribution results in a complete loss of hyperuniformity in both number variance calculations. The inhomogeneities are also reflected in the Fourier- and reciprocal-space images with the emergence of additional but weaker peaks (insets of figure 4(c)).

We performed numerical simulations to clarify the effect of inhomogeneities on hyperuniformity. Figure 5(a) shows the scaled number variance graph of a perfect hexagonal crystal formed by white disks tiling the entire field of view. We introduce voids starting from the corners by removing several disks (figure 5(b)). This action does not affect the number variance calculated through real-space sampling. However, reciprocal-space calculations are immediately affected by this action and show complete loss of hyperuniformity.

The reason for this disagreement is purely geometrical. The observation windows are circular, and the field of view is square (figure 2(b)). These windows cannot sample the corners of the field of view. Consequently, they do not 'see' the voids unless they are extended towards the sampling region (figure 5(c), red square). In contrast, the reciprocal-space calculations immediately detect these voids and loses hyperuniformity since the Fourier calculation effectively projects the real-space onto a torus where the voids overlap and create a large space devoid of white disks.

Figure 5(d) shows a line version of the artificial lattice imperfections directly detected by the reciprocal-space calculations. In contrast, the real-space number variance curve is unaffected by it despite the observation windows having large radii that can sample from the sides of the field of view and 'see' this relatively thin black line (figure 2(b)). However, the line is not sampled by the majority of the windows, and statistically, this imperfection has a negligible effect on the number variance calculations. If we increase the thickness of this line (figure 5(e)) real-space calculations also lose hyperuniformity for large window radii, R > 8. Further increase of the line thickness (figure 5(f)) results in a complete loss of hyperuniformity for real- and reciprocal-space calculations.

The system regains its hyperuniformity at around t = 27 min once the particles are distributed isotopically throughout the system (figure 4(d)). Real- and reciprocal-space calculations both point to the logarithmic, or class II hyperuniformity. Representative number variance curves for t = 27.7 min (frame # 24886) are given in figures 3(b) and (c), and the fitting parameters in table 1. Corresponding Fourier- and reciprocal-space images and the structure factor in appendix figure 18. However, before this happens, video 2 and tables 1 and 2 show an early recovery of hyperuniformity for the reciprocal-space calculations starting at around t ∼ 20 min. Although interparticle distances are clearly non-uniform during 19 < t < 27 min, the power spectrum number variance calculations detect an overall degree of order in the system.

In order to verify this, we performed three numerical simulations where the field of view is divided into two halves. First, both halves are tiled with hexagonal lattices with different interparticle distances (figure 6(a)). Second, one half is tiled with a hexagonal lattice and the other with a disordered stealthy hyperuniform configuration (figure 6(b)). Third, one half is tiled with a hexagonal lattice and the other with a disordered non-hyperuniform configuration (figure 6(c)). Then, the area occupied by the right-hand-side configuration is advanced towards the left-hand-side and eventually fills the entire field of view.

In the first simulation, real- and reciprocal-space calculations point to non-hyperuniform and strongly hyperuniform configurations, respectively (figure 6(a), left). As expected, only the reciprocal-space calculations successfully detected the orderliness of the two halves and confirmed the multi-hyperuniformity of the overall configuration. However, when the areas occupied by the two halves are not equal, in other words if the density variations are large, reciprocal-space calculations also indicate a non-hyperuniform behaviour (figure 6(a), middle). This is similar to the strong inhomogeneity case shown in figure 4(c). When a single hexagonal lattice fully occupied the field of view, both calculations point to a strong, class I hyperuniform state (figure 6(a), right).

The second simulation also confirms these observations with two periodic configurations (figure 6(b)). The third simulation further verifies that if one half of the field of view is filled with a non-periodic, non-hyperuniform configuration, then both calculations point to a non-hyperuniform state at all times (figure 6(c)). In light of these simulations and the power spectrum analyses shown in video 2, we conclude that the particle configurations at 19 < t < 27 min are indeed ordered despite the apparent non-uniform interparticle distances with small density variations.

The conclusion drawn from this sub-section is that real- and reciprocal-space number variance calculations are both necessary and complementary to each other to correctly ascertain the hyperuniformity of a driven dissipative system.

Once the system reaches a statistically isotropic particle distribution, real-space calculations show a brief class I behaviour (at t = 28.5 min frame # 25617), then switch to class II (at t = 32 min frame # 28767) and class III (at t = 34 min frame # 30582) hyperuniform states. Finally, after t = 43.7 min, real-space calculations indicate an irrecoverable non-hyperuniform state as the particles are distancing from each other (tables 1 and 2). Also, we suspected that there is an insufficient number of particles within the field of view to perform reliable real-space hyperuniformity analyses. To test this, we changed the size of our field of view (FOV) for the particle configurations observed at time t = 43.8 min and frame # 39440 (appendix figure 19(a)). The original FOV size is changed from 720 × 720 pixels to a smaller FOV with 680 × 680 pixels and larger FOVs with 800 × 800 and 970 × 970 pixels. The sampling region is kept constant at R_max = 10 μm. Real- and reciprocal-space number variance analyses are presented in table 2, which shows even if the number of particles are increased twice for the FOV size of 970 × 970 pixels, real-space measurements points to the non-hyperuniform configuration.

In addition to the real-space measurements, we also tested the accuracy of our reciprocal-space measurements. We further divided the original FOV to four nonintersecting FOVs with 360 × 360 pixels (appendix figure 19(b)). We measured the reciprocal-space number variance in the quadrants shown in appendix figure 19(b) for R_max > 10 μm, which corresponds to approximately half of the side length. The results provided in table 2 show that the FOV sizes larger than the original FOV do not significantly change the reciprocal-space number variance measurements. However, the FOV sizes much smaller than the original FOV leads to artificial class I or class II hyperuniform states for reciprocal-space measurements. Therefore, we conclude that our original FOV size is sufficient to draw reliable conclusions from the reciprocal-space measurements.

A temporary loss of hyperuniformity for reciprocal-space calculations between 32.4 < t < 32.7 min (frames 29202 to 29451) is due to the detection error for particles at the upper corner of the field of view, which temporarily appear darker with respect to the others. This error can be seen through the dynamic Lindemann parameter map. We note that this frame span is the only such incident throughout video 2. We further note that in figures 4–6, and tables 1 and 2, the density fluctuations scale as ${\sigma }_{N}^{2}\left(R\right)\sim {R}^{x}$ with x > 2, only in the presence of strong inhomogeneities in the system. This is because x = 2 equality only holds for homogenous particle configuration in two-dimensional space.

In contrast to real-space calculations, reciprocal-space calculations continue to show oscillating hyperuniformity between class II and class III states until t ∼ 65 min, and afterwards exhibit only class III hyperuniformity until the end of the video (see figures 3(a) and (b), and tables 1 and 2). This is also reflected in the power spectrum analysis shown in video 2, which indicates a relatively weak periodicity in the particle configuration. To understand this, we analysed the initial frames of the experimental video, right after the sample preparation step and before we turned on the laser for the first time. We repeat this for several experiments. The results for the experiment we report here are shown in the appendix figure 20, where the reciprocal-space image and the structure factor are pointing at a random particle configuration with no periodicity, noticeably different than those seen at the end of video 2.

In our previous study [42], we studied the aggregation mechanism of this particular system. We showed that various interaction mechanisms specific to different particles (CdTe quantum dots and polystyrene particles) and living organisms (E. coli and M. luteus bacterial cells and yeast and human cells) are negligible. In fact, the dominant mechanisms are hydrodynamic and hard-sphere in origin. Therefore, we suspect that the lingering hyperuniformity in this driven dissipative system, long after the absence of an external energy source, may be related to the mutual coupling of the particle positions and velocities with the flow fields.

The hyperuniformity of the system is remarkably robust against random vectorial flows, which can be elastic or nonelastic as well as short or long-lived. Figures 7 and 8 displays time-lapse microscope images during elastic and nonelastic flows, respectively. In elastic flows, the particles do not change their original positions (figure 9(a)), while nonelastic flows sweep clusters of particles and change their original positions (figure 9(b)). New particles then replace these empty sites. Regardless of these differences, the system's hyperuniformity recovers once these flows diminish.

A shrinking bubble boundary forms a short-lived vectorial fluid flow, which can be observed through the particle speed and velocity maps shown in figure 7 and video 2. Despite the apparent inhomogeneities in particle speeds and velocities, the number variance curves are not significantly affected by such flow fields. Similarly, prolonged vectorial flows shown in figure 8 increase the density fluctuations, but then hyperuniformity recovers once they vanish. The flows shown in these particular examples are representative of many similar vectorial flows seen throughout video 2. We believe that the oscillations of the number variance curves observed throughout video 2 are due to perturbations from such short-lived and prolonged vectorial flows.

The robustness of the system's hyperuniformity against various external perturbations directly results from its internal dynamics. As we explained in detail in references [41, 42], the entire system is mutually coupled through intrinsic positive and negative feedback mechanisms. Fluid flows change the particle positions and velocities, and in return, rearrangement of the particles modify the flow field. This is why several lattice imperfections emerge and propagate but do not significantly change the lattice's overall periodicity. Of course, there are limits to such self-correction mechanisms since eachconfiguration has a finite basin of attraction in the phase- and real-space. Sufficiently large spatiotemporal perturbations may well reconfigure the entire system [41]. However, this is not the case in the experiments provided here since the laser was turned off and we have not actively modifying the system dynamics as we previously did in references [41, 42].

In the previous sub-sections, we provided evidence of several hyperuniformity classes through real- and reciprocal-space measurements. Although anomalous suppression of density fluctuations is observed for disordered particle configurations, their structure factor values remain nonzero at infinite wavelengths, $S\left(k\to 0\right){ >}0$. This contradicts the hyperuniformity description established for infinitely large systems, $S\left(k\to 0\right)\sim {k}^{\alpha }$ [2].

Observation of a positive value for structure factor at long wavelengths is not new; it has been reported previously for several experimental systems [6, 26, 53]. It has been argued that the finite size of an experimental system prevents achieving $S\left(k\to 0\right)\sim 0$, even when the system presents clear signatures of hyperuniformity. Dreyfus et al came to the conclusion that reciprocal-space measurements fail to diagnose the hyperuniformity of an experimental disordered jammed soft polydisperse particles, where real-space measurements are deemed more reliable [6]. Berthier et al confirmed the hyperuniformity of a jammed packing of size-disperse spheres by looking at the vanishing isothermal compressibility instead of the static structure factor [26].

In addition to the finite size effects, the effects of finite-resolution of particle positions, experimental and measurement noise, imperfections in optical imaging, image size restrictions, and size polydispersity are also pointed out as difficulties to diagnose the hyperuniformity of an experimental system directly from structure factor analysis [6, 26]. Admittedly, all these effects limiting the potential of an experimental system accessing to small k values are also present, to some extent, in our experimental system, except size polydispersity. However, unlike in previously reported cases, our real- and reciprocal-space measurements mostly agree on this particular system's hyperuniformity. Earlier studies measure hyperuniformity directly from the power spectrum. In addition to this, here we also provide the number variance of the power spectrum (4). This fact, together with our results discussed in the previous sections led us to conclude that both measurements have their merits. Real- and reciprocal-space number variance measurements complement each other and should be used together to assess the hyperuniformity of a given experimental system correctly.

We believe that the physical reason for the agreement of our dual approach on the system's hyperuniformity is the mutual coupling of the entire system, as pointed out in the previous sub-section. Further, the positive structure factor value for small k values is likely due to the Brownian motion of the particles and structural defects produced by their random collisions in our experiment.

We tested both of these hypotheses with numerical simulations. The results are presented in figure 10. A dynamic Lindemann parameter map of the experimental frame at t = 28 min (frame #25184), the corresponding scaled number variance curves, and the structure factor are shown in figure 10(a). The real- and reciprocal-space calculations point to a class I and class II hyperuniform states, respectively with $S\left(k\to 0\right)\sim 0.105$. In figure 10(b), we show the same analyses for a simulated stealthy hyperuniform structure exhibiting class I hyperuniformity with $S\left(k\to 0\right)=0$. Next, we introduce Brownian noise to each particle in the stealthy configuration (figure 10(c)), for which the number variance curves indicate a class II hyperuniform state with $S\left(k\to 0\right)\sim 0.09$. Noticeably, the peaks at k > 5 vanish, pointing to the presence of only short-range interparticle correlations. We then introduce several structural defects without the Brownian noise (figure 10(d)), which, again, results in a class II hyperuniform state with $S\left(k\to 0\right)\sim 0.128$. Lastly, we include both Brownian motion and structural defects (figure 10(e)). The system exhibits a class II hyperuniformity with $S\left(k\to 0\right)\sim 0.099$. The corresponding structure factor graph similarly indicates short-range interparticle correlations as in figure 10(c).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We believe the Brownian motion of the particles and the structural defects produced by their collisions is sufficient to explain the positive structure factor for small k that we observe experimentally. We also showed that the structure factor of the experiment showing oscillating peaks indicating long-range particle correlations, which are likely preserved thanks to the nonlocal interactions ensured by the flow fields.