{'fraction_non_alphanumeric': 0.02883645240032547, 'fraction_numerical': 0.00836452400325468, 'mean_word_length': 4.319598337950138, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Kinetically-constrained freezing transition in a dipole-conserving system

We study a stochastic lattice gas of particles in one dimension with strictly finite-range interactions that respect the fracton-like conservation laws of total charge and dipole moment. As the charge density is varied, the connectivity of the system's charge configurations under the dynamics changes qualitatively. We find two distinct phases: Near half filling the system thermalizes subdiffusively, with almost all configurations belonging to a single dynamically connected sector. As the charge density is tuned away from half filling there is a phase transition to a frozen phase where locally active finite bubbles cannot exchange particles and the system fails to thermalize. The two phases exemplify what has recently been referred to as weak and strong Hilbert space fragmentation [1, 2], respectively. We study the static and dynamic scaling properties of this weak-to-strong fragmentation phase transition in a kinetically-constrained classical Markov circuit model, obtaining some conjectured exact critical exponents.

Refs. [1,2,[56][57][58] showed that the combination of strictly finite-range interactions and the fracton-like constraints of charge and dipole conservation results in a fragmentation of Hilbert space into exponentially many (in volume) dynamically disconnected sectors, herein referred to as Krylov sectors. This means that a graph, where nodes represent charge configurations and edges represent local dipole-conserving transitions, consists of exponentially many disconnected components within each of the polynomially many symmetry sectors of configurations with a common charge and dipole moment.
Hilbert space fragmentation comes in two distinct types: strong and weak [1,2]. In systems that are strongly fragmented typical initial charge configurations fail to thermalize because they are constrained to explore only a vanishingly small fraction of the states with the same charge and dipole moment. On the other hand, in the case of weak fragmentation the configuration space still shatters into exponentially many disconnected Krylov sectors, but the fraction of states that belong to the largest sector approaches one in the limit of large systems, and therefore typical initial states are in this largest Krylov sector and thus can potentially thermalize.
In this paper we study a charge and dipole-conserving model that is weakly fragmented at charge densities around half filling, for which typical initial states do thermalize. However, as the total charge in the system is varied, a critical point is encountered beyond which thermalization cannot occur due to strong fragmentation of the set of all charge configurations with that total charge. Since the fragmentation in this case is dictated by fundamentally classical constraints-because we are working in a basis where fractonic degrees of freedom are assumed to exist without having to emerge from an underlying quantum model-we study a kinetically-constrained [59,60] classical Markov circuit model. We study the static and dynamic properties of the phase transition between the weakly and strongly fragmented phases by constructing and supporting a simplified theoretical model that we believe captures many of the key features of the critical point, as well as by simulating the full dynamics. Our simplified model gives conjectured exact critical exponents of ν = 2 for the correlation length, and β = 1 for the density of frozen sites which serves as an "order parameter" for the frozen phase. The dynamics in the thermalizing, weakly-fragmented phase away from the phase transition is demonstrated to be subdiffusive, with the transport time growing as the fourth power of the length. At the phase transition, the dynamics is slower, as expected, with a dynamical critical exponent z that can shuffle the local charge configuration in any way that conserves the total charge and dipole moment. Hollow black circles represent the charge configuration before the action of the gates, and pink circles represent it afterwards. Each site can have 0, 1, or 2 particles. The black arrows represent possible actions of the gates; the only transitions that are allowed in this model are "pair hopping" moves where one (or two) particles hop to the left and one (or two) hop to the right. This model has a detailed balance so that in the longtime steady state all moves are equally likely to occur as the reverse of the same move.
appears to be near 7, but this can only be explored numerically over a modest range of length scales, so the true asymptotic z might be larger than this.

II. MODEL
We study a model of indistinguishable particles hopping on a one-dimensional lattice of L sites. Interactions are strictly k-local and each site i can host 0 ≤ n i ≤ n max particles. Equivalently, one can consider classical spin models where each site hosts a spin S = (n max − 1)/2 degree of freedom, and particle configurations map to spin configurations in the z basis. Both the total charge (particle number) N 0 = i n i and dipole moment N 1 = i (x i − x 0 )n i are conserved by the dynamics, where x 0 can be chosen to be the midpoint of the system, and the sites are at integer positions 1 ≤ x i ≤ L. The dynamics are given by a kinetically-constrained classical Markov circuit. The circuit is made up of layers of k-local Markov gates that stochastically map any local charge configuration to any other with the same charge and dipole moment on those k consecutive lattice sites. A time step consists of one tightly-packed layer of gates with a randomly chosen spatial shift of up to k−1 sites. Most of our results are for systems with open boundaries. Gates that extend past either edge of the system must act only on less than k sites. Fig. 1 depicts one step of the time evolution. We use gates that assign an equal probability to all possible transitions, i.e. the "infinite temperature heatbath" case. The average of any quantity in the steady state is then given by averaging over all charge configurations in the Krylov sector of the initial state. Doing this average exactly is only possible for small systems in practice, so we statistically sample the distribution for larger systems for which we can simulate the steady state, while for even larger systems near the phase transition we can only numerically access nonequilibrium dynamics. For the rest of this paper we focus on the specific case of n max = 2 and k = 4, which connects the charge configuration space in exactly the same way as the 4-local spin-1 quantum models of Refs. [1,2], and on systems that have approximately zero total dipole (N 1 = 0). Thus the phase diagram we study is as a function of the average charge densityn = N 0 /L, which is the parameter that tunes the system from weakly to strongly fragmented. Systems with n max = 2 and k < 4 do not have a weakly fragmented phase at all, so k = 4 is the minimal generic case for n max = 2. If one considers systems with a substantially nonzero dipole moment, the equilibrium charge density is spatially nonuniform and the freezing transition occurs locally where the local charge density reaches its critical value. For these systems with zero total dipole, we say the system thermalizes if in the limit of a large system and long times for almost all initial states the time-average of the local density isn at all sites that are far from the ends of the system.

A. Exact enumeration
We first study small systems (L ≤ 18) for which we can eliminate finite-time effects, but not finite-size effects, by enumerating all possible configurations and constructing all of the Krylov sectors exactly. This is similar to the approach of Refs. [1,2]. The authors of those works defined a quantity to diagnose strong fragmentation that we herein call D max /D sum ; it is the ratio of the number of configurations in the largest Krylov sector (D max = max j D j ) to the total number of configurations in all of the Krylov sectors combined (D sum = j D j ). Since charge and dipole are conserved, the index j runs over all Krylov sectors corresponding to a certain total charge N 0 and dipole N 1 of interest, and D j is the number of charge configurations in sector j. To make contact with previous studies we plot this quantity for a range of small system sizes and charge densities in Fig. 2(a). This shows an onset of strong fragmentation as the charge density is tuned away from half filling; however finite-size effects are strong and we cannot argue from these data alone that there is a sharp transition in the thermodynamic limit.
In an attempt to improve on this measure of fragmentation in small systems we define the entropy of fragmentation where p j = D j /D sum and j max is the total number of Krylov sectors being summed over. We do this to minimize the effects of D sum becoming small asn → n max in very small systems; these effects cause the quantity in Fig. 2(a) to curve downwards asn → 2. A weakly fragmented system is characterized by lim L→∞ s frag = 0, and a strongly fragmented one will have s frag > 0. The entropy of fragmentation is shown in Fig. 2 Both of the diagnostics of strong fragmentation discussed so far rely on constructing the Krylov sectors exactly, which is prohibitively expensive for all but very small systems. Furthermore the charge density at these small system sizes is not an approximately continuous tuning parameter and this adds difficulty to studying the phase transition. We therefore propose a different diagnostic that can be computed from local measurements and thus can be estimated in larger systems with some finite-time effects without constructing the full Krylov sectors. It is inspired by studies of kinetically-constrained classical models of glasses [61]: We propose that the frozen, strongly fragmented phase can be characterized by a nonzero fraction of frozen sites in a typical history, where we say that site i is frozen at time t if n i has remained unchanged for all times ≤ t. As time evolves frozen sites can become active, but with this definition active sites never become frozen. This is the same definition of active and frozen sites used in the discussion of "shielding regions" and the construction of nonthermal states in Refs. [1,2,56]. In those works the focus was on systems at half filling with various ranges of interactions k, and on constructing the partially frozen states even in weakly fragmented systems where these are vanishingly rare. To relate the notion of frozen sites back to Krylov sectors and fragmentation: for a given initial state a site will remain frozen for all times if for every configuration in its Krylov sector that site has the same charge n i . The fraction of frozen sites ρ F at infinite times, averaged over all initial states, serves as an order parameter heralding the onset of the strongly fragmented phase. Much of our analysis is focused on this notion of frozen and active regions. In Fig. 2(c) we show the infinite-time ρ F averaged over all states with charge densityn and total dipole N 1 = 0. This shows finitesize indications of the phase transition between weakly and strongly fragmented phases near the nominal critical densityn c = 1.5, where ρ F becomes nonzero.

B. Active vs. frozen regions
In the vicinity of the critical density large systems typically contain regions that are initially frozen and regions that are initially active. The regions that are initially frozen are typically contiguous blocks of sites with a local charge density near n max = 2, and the initially active regions are typically contiguous blocks of sites where the charge density is closer to half filling and below the critical density. For example, by identifying the initially frozen sites as the sites that cannot evolve to any other value of charge during the first step of time evolution, no matter what spatial shift that first layer of the Markov circuit has, we find that systems with L = 10 6 and n = 1.5 contain contiguous blocks of frozen (F ) sites with mean charge density n F ∼ = 1.95 and active (A) blocks with n A ∼ = 1.28, where the averages are taken over equally weighted blocks. As time goes on, a frozen site can become active if a neighboring active block is a good enough charge bath to move some of its charge, so the active blocks can thus grow and merge. As the active blocks grow and occupy a larger fraction of the system their charge density increases. The growth stops when the charge density of the active block approaches the critical value and the active block thus stops being able to move the charges on the neighboring frozen sites.

C. Single active blocks
In order to understand this process we study the idealized scenario of a half-filled initial block of L 0 contiguous sites, with initial charges n i = 1 on all of those L 0 sites, embedded in an otherwise fully-filled system of infinite length with initial charges n i = 2 on all other sites. This is in some ways similar to the idealized scenarios used to study the effects of thermal inclusions in an MBL system and their relevance to the MBL transition [24, 26-28, 62, 63]. As the system evolves we keep track of which sites have become active and the corresponding charge density n A of the growing active bubble. When the length of the block of active sites has grown to L A , then the charge density of the active block is n A = 2 − (L 0 /L A ). In Fig. 3(a) we show the resulting time evolution of the sample-averaged n A for active bubbles of initial size L 0 = 12 to 32. In each such sample with even L 0 ≤ 48 the active bubble converges to exactly twice its initial length and a charge density of n A = 1.5 ∼ =nc, where it finally ceases to be a good enough charge bath to move charges from the neighboring frozen sites. We have been able to show by an explicit iterative construction that even length bubbles can grow to twice their length, thus showing thatn c ≥ 1.50; we suspect that the exact value of the critical density does saturate this bound, but we have not yet been able to prove that.
We have also studied initial blocks of (L 0 + 1) contiguous sites where L 0 of the sites in the block are singly occupied, while one of the sites within the block and all sites outside of the block are doubly occupied. Such active bubbles all appear to stop their growth with length 2L 0 contiguous active sites, with the sole exceptions of all cases where the center of mass of the active bubble is precisely on a site. These latter cases instead stop with length (2L 0 − 1). We have numerically examined many bubbles with even and odd L 0 ≤ 48 and various locations of their centers of mass, not finding any exceptions to this behavior.
We denote the time at which the active bubble first reaches its final length by τ and extract a dynamical scaling exponent z ∼ = 7.0(5), such that τ ∝ L z 0 , by collapsing the distribution of τ over samples for each system size, as shown in Fig. 3(b). However the range of L 0 over which we were able to compute τ is quite small (on a log scale), and statistically significant drift towards higher z is present (see Fig. 3(c)), thus it is possible that z drifts upwards indefinitely with increasing L 0 , or it may converge to a finite value. We therefore leave our estimate of z ∼ = 7 as a lower bound. We also performed the same "ideal active bubble" experiments but with active blocks that were initialized in a random configuration at half filling embedded in a fully-filled system. In that case the density of the active blocks of almost all samples still converges to n A = 1.5, and the fraction of samples for which this is not true decreases with L 0 . We therefore expect that the ideal scenario studied above is representative of large active bubbles embedded in typical large systems. In a typical system there will be many such active blocks, and they will expand and merge as time goes on. We propose that these blocks expand until they self-tune to the critical charge density where they are no longer good charge baths, and the critical point is where the system is just able to be completely covered by the largest, slowest active block at the latest times. Just above the critical density, a finite density of frozen sites will persist indefinitely and block charge transport, preventing thermalization.

D. Approximate model for the transition
By considering isolated active bubbles embedded in otherwise frozen systems we have developed the idea that an active bubble typically unfreezes the immediately adjacent frozen sites and, in doing so, increases its average charge density and length. This process continues until the critical charge density is reached, at which point the growth stops because the active bubble is no longer an effective charge bath. We therefore propose an ap- proximate scheme for partitioning a system withn >n c into frozen and active blocks at time t = ∞: Given a system at t = 0 we first compute its initially active (A) and frozen (F ) blocks of sites by considering which sites could evolve to a different charge during the first time step, regardless of the spatial shift of that first layer of gates in the Markov circuit. We then loop over all of the active blocks, allowing each block to expand by incorporating the sites and the charges in neighboring frozen blocks. This expansion proceeds in either or both directions and is continued as long as the active block's average charge density satisfies n A <n c and it does not overlap with another active block. After allowing all of the active blocks to grow like this, any two active blocks that have come into direct contact are merged and considered as one block from then on. We iterate this process of growing and merging until it has converged, i.e. all active blocks reach density n A =n c and are thus unable to grow farther. This approximate model of the steady state must be given the critical density as an input parameter, so we choose to usen c = 3/2, which is our current best (and perhaps exact) estimate of the critical point. For n <n c this model will instead grow one active block that encompasses the entire system. Although we will provide analytic arguments for all of the following results, it is convenient to first present the numerical data from simulations of this approximate model of the steady state. At each of a series of charge densitiesn approachingn c from above we simulate this process on 10 2 samples of large systems (L = 10 6 ), initialized randomly from all configurations with the correct amount of total charge, and store the distribution of the lengths L A of the final active blocks; the final frozen blocks remain of length of order one forn near but abovē n c . The resulting data are shown in Fig. 4(a). The distributions, which we denote by pn(L A ), approach a power law ∼ L −α A with α = 3/2 asn approachesn c from above. We define a correlation length ξ = L 2 A / L A , which we compute from the data at each value ofn, and we plot this in Fig. 4(b). These averages are over all active blocks, with each block equally weighted independent of its length. Defined in this way, ξ is a characteristic active block length where the distribution changes from a power law for L A ξ to something steeper for L A ξ. From this we see that ξ diverges at the critical density as ξ ∝ (n−n c ) −ν , with ν = 2 being consistent with the data. We assume the scaling ansatz pn( where f is a universal function, to collapse the data of Fig. 4(a) and we show the scaled distributions in Fig. 4(c). With this scaling, L p A ∼ ξ (p+1−α) for any power p > α − 1, which is why we defined ξ the way we did, rather than as L A , for example.
To summarize: This approximation yields the predictions that the critical point is characterized by a powerlaw distribution of the lengths of active blocks with an exponent of α = 3/2, and by a correlation length ξ which diverges with the exponent ν = 2 as the critical point is approached from the strongly fragmented phase.
In order to argue analytically for the critical powerlaw distribution pn c (L A ) ∝ L −3/2 A we consider a critical (n =n c ) system with an embedded active block whose total charge is initially ∆N 0,A below what it would be if the block were at the critical density. The block will begin to grow and each time its length increases by ∆L A it gains an amount of charge ∆L Anc + δ, where δ is a random variable with mean zero. Therefore ∆N 0,A , the total charge deficit of the block relative to the critical total charge, does an unbiased random walk in one dimension with an "absorbing wall" at ∆N 0,A = 0 where the active block becomes critical and can no longer grow. The analogous real-space problem is that of a random walk x(t), where x is position and t is time, that begins at a negative position x(t 0 ) = −x 0 and we want to know: What is the distribution of final times t f at which the walker reaches x = 0 for the first time? This can be solved in the continuum by the method of images, with the result being that the distribution of final times goes like p(t f ) ∼ t −3/2 f , which is analogous to the observed result that pn c (L A ) ∼ L −α A with α = 3/2. We can also argue for the result ν = 2: The frozen blocks of the system have an average charge density n F that is always well aboven c and close to n max , and we have argued that large active blocks will converge to the critical charge density n A =n c at late times. Therefore we consider an approximation where all frozen blocks have n F = n max = 2 and all active blocks have n A =n c as t → ∞. In that case the shifted charge density can be expressed asn −n c = (n max −n c ) L F /( L A + L F ), and since L F does not diverge at the transition this implies L A ∝ (n −n c ) −1 as the critical point is approached within the frozen phase. If the distribution of L A takes the form pn ∝ L −3/2 A f (L A /ξ), then the identity L A = pn( )d and dimensional analysis implies L A ∝ √ ξ and therefore ξ ∝ (n −n c ) −2 , i.e. ν = 2. We can also obtain the critical exponent β of the order parameter ρ F ∼ (n −n c ) β (the fraction of frozen sites) via similar reasoning: Again because L F does not diverge at the transition and L A diverges as (n −n c ) −1 , the identity ρ F = L F /( L F + L A ) implies β = 1.

E. Simulation results
At this point our understanding has been built up from idealized initial conditions and simplified models. We now go back to the full dynamics of large systems in an effort to validate as much of the picture we developed above as possible. We randomly initialize systems with definite charge densityn of size L = 10 4 and run their full dynamics up to t = 10 8 layers of the Markov circuit. This is still very early times for such large systems (see Fig. 3), so we are targeting the infinite system, finite time limit, contrary to the infinite time, small L results shown in Fig. 2. In Fig. 5 we show the time evolution of the fraction of frozen sites ρ F and the charge density of the active fraction of the system n A . Fig. 5(a) is consistent with the result of Fig. 2(c) for small systems and it appears that ρ F ∝ (n −n c ) β with β = 1 andn c = 1.50 is consistent in the strongly fragmented phase, as our simplified theory and previous numerical results suggested. The slight deviation of the late-time data from the simplified theoretical prediction (dashed line) deep in the frozen phase is understood: In typical systems the frozen blocks contain some sites that do not have charge n i = n max , thus ρ F → 1 slightly beforē n → 2, and indeed this was observed in small systems as well (see Fig. 2). This effect is significant deep in the frozen phase because that is where the frozen blocks are largest and can contain a small number of isolated frozen sites that have n i < n max .
In the second plot, Fig. 5(b), we show that the charge density of the combined active blocks of the system is converging to n A =n from below in the active phase (low density active blocks activate high density frozen blocks until the system is entirely active) and in the frozen phase near the critical point there is an extended regime where the active portion of the system self-tunes to the critical density n A →n c , as we have suggested. Deep in the frozen phase almost all of the active bubbles are very small and many of them have a two-state Krylov sector (...2220222... and ...2212122...) that grows only to density n A = 4/3, so the average n A becomes noticeably less thann c (dashed line), as is visible in Fig. 5(b). Finally in Fig. 5(c),(d) we show how ρ F and n A are converging to their final values at eachn in order to demonstrate that n c = 1.5 is a reasonable estimate of the critical density for typical large systems. Since we are only able to access rather early time dynamics, we have not been able to test the prediction that the critical distribution of the lengths of active blocks goes like pn c (L A ) ∝ L −3/2 A . This may be because we are not able to build up large active blocks from small ones during the full simulation of large typical systems during the accessible times. Thus we leave the results presented in Fig. 4 as predictions of our simplified theory, not yet tested by simulations of the exact dynamics.

F. Subdiffusion in the thermalizing phase
Finally, much of our analysis has focused on the properties of the strongly fragmented phase, and we have viewed the transition as the breakdown of strong fragmentation as the charge density of the system is tuned closer to half filling and active blocks are able to thermalize their surrounding frozen neighbors. Before concluding, here we establish a basic property of the dynamics Autocorrelation of charge density fluctuations at half filling. The dotted line represents the scaling C(t) ∝ t −1/4 . Data for C(t) (black crosses) was generated by simulating the dynamics of 10 3 samples of systems with L = 10 3 sites and periodic boundaries out to time t = 10 6 . Averaging was done over all sites, samples, and windows of time surrounding each data point shown above.
of charge in the thermalizing, weakly fragmented phase. In Fig. 6 we show the averaged autocorrelation function of the charge fluctuations C(t) = (n i (t) −n)(n i (0) −n) for systems with L = 10 3 ,n = 1, and periodic boundaries. These are the only simulations for which we employ periodic boundaries because it allows for averaging over all sites. The average was taken over all 10 3 sites of the system, 10 3 samples, and over the intervals of time shown in Fig. 6. We find that the charge dynamics of the weakly fragmented phase is subdiffusive with a characteristic space-time (x-t) scaling x ∝ t 1/4 , and we believe this dynamical scaling is characteristic of dipole-conserving systems.
An argument for the subdiffusive x ∝ t 1/4 scaling goes as follows: We work with a one-component hydrodynamics of the charge density n(x, t), which is now a continuous-valued field in the coarse-grained system. Since the dipole moment is conserved, the fundamental current is not one that displaces charge from position x to x + dx, instead it is one that displaces equal amounts of charge from x to x + dx and x − dx, and thus preserves the dipole moment. We call the density of these processes the pair hopping density p(x, t), and take the convention that a positive p(x, t) depletes the charge at position x. Since these processes are isotropic, a gradient of them can generate a local charge current j, i.e. j ∝ −∇p. Assuming the system is entropy-driven and that the entropy density at position x only depends on the density at that position, the change in entropy due to local pair hopping is ∝ −∇ 2 n when expanding around a uniform-density equilibrium, and thus charge currents are driven by the "Fick's Law" j ∝ ∇ 3 n. Combining this with charge conservationṅ + ∇j = 0 yields the desired resultṅ ∝ −∇ 4 n, which implies an x ∝ t 1/4 scaling for the relaxation of charge density perturbations.
This subdiffusive scaling is also in agreement with a more complete theoretical framework for the hydrodynamics of thermalizing systems of fractons that was recently developed [64], and it was also recently observed in a cold atom quantum simulation of a "tilted" Fermi Hubbard model [65] with an emergent dipole moment conservation on lengthscales larger than the lattice spacing. Such a platform could provide an experimental testing ground for studying weak and strong Hilbert space fragmentation in quantum systems, and for investigating phase transitions between the two like the one studied in this paper.

IV. SUMMARY AND OUTLOOK
We studied the freezing phase transition encountered when tuning the average charge density of a onedimensional system that conserves its total charge and dipole moment. In the thermalizing phase the system is weakly fragmented in the sense that almost all charge configurations are connected to each other by the stochastic dynamics of our model, and thus the system eventually reaches a global equilibrium from typical initial states. In the frozen phase the set of global charge configurations is strongly fragmented, so the set of charge configurations shatters into exponentially many dynamically disconnected sectors, with no single sector being dominant over all others, and we explained how this manifests itself in terms of the growth and isolation of locally weakly fragmented active bubbles. These bubbles act as charge baths and help to thermalize the initially frozen regions of the system, but in the strongly fragmented phase active bubbles eventually stop growing and they remain isolated from each other for all time. We studied the dynamic scaling of these active regions and found that they self-tune to the critical charge density near the critical point. Based on this understanding we developed a solvable simplified model that yields several predictions about the static properties of the critical point. We also showed that in the weakly fragmented phase the charge dynamics is subdiffusive, and discussed the prospects for studying these dipole-conserving systems experimentally. This work has focused on the classical aspects of Hilbert space fragmentation in dipole-conserving systems. Two related question for future research are then What are the distinctly quantum aspects of systems that conserve dipole moment? and How is the physics of dipole-conserving systems related to the glassy phenomenology of other classical kinetically-constrained models?
Note added: While finalizing this manuscript we became aware of related work by J. Feldmeier, P. Sala, G. de Tomasi, F. Pollmann, and M. Knap. Where our results overlap they agree.