{'fraction_non_alphanumeric': 0.044937342597930634, 'fraction_numerical': 0.019783354067421186, 'mean_word_length': 4.209927284223838, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 6, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Critical fragmentation properties of random drilling: How many random holes need to be drilled to collapse a wooden cube?

A solid wooden cube fragments into pieces as we sequentially drill holes through it randomly. This seemingly straightforward observation encompasses deep and nontrivial geometrical and probabilistic behavior that is discussed here. Combining numerical simulations and rigorous results, we find off-critical scale-free behavior and a continuous transition at a critical density of holes that significantly differs from classical percolation.

A solid wooden cube fragments into pieces as we sequentially drill holes through it randomly. This seemingly straightforward observation encompasses deep and nontrivial geometrical and probabilistic behavior that is discussed here. Combining numerical simulations and rigorous results, we find off-critical scale-free behavior and a continuous transition at a critical density of holes that significantly differs from classical percolation. The connectivity of a solid block of material strongly depends on the density of defects. To systematically study this dependence one must find an experimental way to create defects inside the solid. For example, in 2D one can simply punch holes in a sheet and measure the physical properties of the remaining material. But in 3D, inducing localized defects is not simple. One conventional solution consists in perforating the material by drilling holes or laser ablation from the surface [1,2].
In a table-top experiment, we start with a solid cube of wood and plot on each face a square-lattice mesh of L by L cells. Initially, the cube has no holes. Sequentially, for each one of three perpendicular faces, we randomly choose one square-cell and drill a hole having a radius of 1/ √ 2 cell lengths to the other side of the cube. We repeat this process iteratively until the entire structure collapses into small pieces and the bottom and top part of the cube are no longer connected. The first row of Fig. 1(a)-(d) shows the result of the drilling process of a real cube with edge length 6 cm (manufactured from 2 cm thick plates of medium-density fiberboard (MDF)), where holes were drilled with a diameter of 1 cm. As the drilling proceeds, pieces get disconnected and eventually the entire structure collapses.
Numerically, we start with a three-dimensional cubic lattice of L 3 sites and fix three perpendicular faces. A fraction 1 − p square-cells on each face is randomly selected and all sites along the line perpendicular to that face are removed (second row of Fig. 1). Thirty years ago, Y. Kantor [3] numerically studied this model on lattices of up to 10 6 sites and concluded that the critical fragmentation properties of this model are in the same universality class as random percolation [4][5][6]. Here, we combine rigorous results and large-scale numerical simulations, considering lattices three orders of magnitude larger in size, to show that this is not the case. Removing entire rows at once induces strong long-range directional correlations and the critical behavior departs from random percolation. Also remarkably, while in random fragmentation power-law scaling is solely observed around the critical threshold, here we find it in an entire off-critical region. These findings suggest that longrange directional correlations lead to a rich spectrum of critical phenomena which need to be understood. Possible implications for other complex percolation models are discussed in the conclusions.
Threshold. The average total number of drilled holes is 3(1 − p)L 2 and the asymptotic probability that a site in the bulk is not removed is p 3 . We first measure the threshold p c at which the cube collapses for different lattice sizes, up to L = 1024, using different estimators of the transition point, as discussed in the Supplemental Material [7]. Extrapolating the data to the limit L → ∞, gives p c = 0.6339 ± 0.0005, consistent with the value estimated by Kantor using Monte Carlo renormalization group techniques [3] (see Supplemental Material [7]). This threshold is larger than the two-dimensional squarelattice percolation threshold (p 2D ) [8,9] and smaller than the cubic root of the one for the three-dimensional simple cubic lattice [10]. Static exponents. We consider the fraction P ∞ of sites in the largest cluster of connected sites (see Supplemental Material for more data of P ∞ (p) and of other observables [7]). P ∞ is the standard order parameter in percolation identifying the transition from a disconnected to a globally connected state. For the drilling model, the situation will turn out to be more complicated. Figure 2(a) shows a double-logarithmic plot of the order parameter, rescaled by a power of the lattice size P ∞ L β/ν as function of the distance to the transition |p − p c |L 1/ν . Based on finite-size scaling analysis [5], we find that the critical exponent of the order parameter is β = 0.52 ± 0.04, and the inverse of the correlation length exponent is 1/ν = 0.92 ± 0.01 (see Supplemental Material). We note that both β and ν are different from the corresponding values for 2D and 3D classical percolation. However, somehow surprisingly, the exponent ratio β/ν = 0.50 ± 0.04 is within error bars the same as for 3D percolation. Thus, while the fractal dimension of the largest cluster (given by d f = d − β/ν) is consistent with the one for 3D percolation, the larger value of β (compared to 3D percolation) implies that the transition from the connected to the disconnected state is less abrupt (see also Supplemental Material). We consider next the behavior of the second moment of the cluster size distribution M 2 , excluding the contribution of the (color online) (a) Double-logarithmic plot of the rescaled order parameter P∞L β/ν as function of the scaling variable (p − pc)L 1/ν for different lattice sizes L. The linear part has a slope of β = 0.52 ± 0.04, consistent with β/ν being the same as for three-dimensional percolation, but with a different exponent 1/ν = 0.92. The present value of β is different from the two-dimensional one β = 5/36 ≈ 0.139 [5,11] and the one in three dimensions, β ≈ 0.417 [12]. (b) Double-logarithmic plot of the rescaled second moment M 2 L −γ/ν , with γ/ν = 2.0452, as function of the scaling variable |p − pc|L 1/ν , with 1/ν = 0.92, for different lattice sizes L. The solid black line is a guide to the eye with slope −2.3. The two sets of data points correspond to the sub-critical (p < pc) and super-critical (p > pc) regions.
largest cluster. As shown in Fig. 2(b), the finite-size scaling analysis gives that the susceptibility critical exponent is γ = 2.3 ± 0.1. We note that our results for the static critical exponents are within error bars consistent with the scaling relation 2β + γ = dν (for d = 3).
Dynamical exponents. The transport properties of the largest cluster at the critical threshold, p c , are intimately related to dynamical critical exponents and they can be measured by quantifying three sets of sites in the largest cluster [13]. First, we consider the so-called red sites.
A site is considered a red site if its removal would lead to the collapse of the largest cluster [14]. The red sites form a fractal set of fractal dimension d RS = 0.92 ± 0.07 (see Fig. 3), which is compatible with the inverse of the correlation length exponent ν that we obtained from the finite-size scaling analysis in Fig. 2, d RS = 1/ν, as predicted by Coniglio for classical percolation [15]. However the value of d RS = 1/ν for the drilling transition is very different from the classical 3D percolation result 1/ν = 1.1437 ± 0.0006 [12]. Figure 3 also shows that the shortest path connecting the top and bottom sides of the largest cluster is a fractal of fractal dimension d SP = 1.30 ± 0.05. Finally, the backbone of the largest cluster between its bottom and top ends is defined as the set of sites that would carry current if a potential difference is applied between the cluster ends (also known as bi-connected component). The backbone fractal dimension is determined as d BB = 2.12 ± 0.08, which is larger than in classical 3D percolation, where d BB = 1.875 ± 0.003 [13,16,17]. Qualitatively, an increase in the backbone fractal dimension is compatible with a simultaneous decrease in the shortest path fractal dimension, since both correspond to a more compact backbone, similar to what is observed in long-range correlated percolation [18,19]. Thus, although the fractal dimension of the largest cluster is similar in both classical percolation and drilling, the internal structure of the largest cluster is significantly different. This implies that transport and mechanical properties of the largest cluster follow a different scaling.
Cluster shape. Given the highly directional nature of the drilling process, we analyze the symmetry of the different clusters. In particular, we consider them as rigid bodies, consisting of occupied sites at fixed relative positions, and look at the eigenvalues and eigenvectors of their inertia tensors [22,23]. The numerical results show that, when compared to classical percolation clusters, the drilling transition clusters are more anisotropic, their orientations being mainly aligned along the direction of the cube edges (see Supplemental Material for quantitative details [7]).
We now give a rigorous argument for the existence of asymmetric clusters in drilling percolation. Fix some p ∈ (p 2D , p c ), where p 2D ≈ 0.5927 < p c is the critical threshold for 2D site percolation. Consider a lattice size L X × L Y × L Z with L X = L Y = L and L Z = e L and take a square domain A in its base with side length k = c 0 log(L) where c 0 is a positive constant that is smaller than − [log (p(1 − p))] −1 . Say that the event S(A) occurred if, along the z-direction, no point inside A is drilled but all points on its boundary are. For large L, this event happens with probability at least L c0 log (p(1−p)) . Consider also two rectangles R x and R y in the (x, z) and (y, z)-planes, respectively, aligned with A. These rectangles have base length k and height exp{c 1 k}, where c 1 is an arbitrary positive constant smaller than (color online) Number of sites in the backbone of the spanning cluster, length of its shortest path (chemical distance), and red sites in the backbone, measured at p = pc = 0.6339, as function of the lattice size L. Considering the local slopes of the data, we obtain the following fractal dimensions: dBB = 2.12 ± 0.08, dSP = 1.30 ± 0.05, and dRS = 0.92 ± 0.07. The fractal dimension of the red sites is, within error bars, compatible with the value 1/ν ≈ 0.915, found from the finite size scaling behavior of the pc estimators (see Fig. S1), and the relation dRS = 1/ν [15]. Here, dBB is larger than in classical three-dimensional percolation, where dBB = 1.875 ± 0.003 [13,16,17], while dSP is smaller than the classical value dSP = 1.3756 ± 0.0006 [20]. The solid lines are guides to the eye. To extract the fractal dimensions, we analyzed the local slopes as proposed in Ref. [21]. Results are averages over at least 5 × 10 3 samples. the correlation length for two-dimensional percolation with parameter p. The event R x (respectively R y ) indicates the existence of a path crossing R x (respectively R y ) from bottom to top that has not been drilled in the y (respectively x) direction. By our choice of c 1 , R x and R y have positive probability (uniformly over k), i.e. there exists a δ > 0 such that P (R x ) ≥ δ. In addition, if S(A) occurs, then there exists a cluster spanning A × [0, e c1k ] from bottom to top and whose projection into the (x, y)plane does not extend beyond A. Thus the probability of finding a cluster of radius k and height e c1k is bounded from below by the probability that there exists a square A along the diagonal x = y for which S(A) ∩ R X ∩ R Y occurs, which is greater than which converges to unity as L increases. This shows that one expects to have clusters extremely aligned along the z axis, as numerically observed (see for example Fig. S15 of the Supplemental Material [7]). In fact, the same argument can be straightforwardly extended to explain the alignment along the x and y directions, as also observed. Spanning probability. To understand the properties of the drilling transition in terms of global connectivity, we consider the spanning probability Π(p), defined as the probability to have at least one cluster including sites from the top and bottom of the lattice, at a given value of the control parameter p. Figure 4 shows the spanning probability below the threshold Π(p = 0.63), for different lattice aspect ratios. The lattice size is L X × L Y × L Z with L X = L Y and L Z = rL X and the spanning probability is measured in the z-direction. At the drilling transition, p = p c , the spanning probability approaches a constant for large lattice sizes (see Supplemental Material, Fig. S9), similar to what is observed for classical percolation [24,25]. By contrast, for values of p between p 2D and p * , the numerical results suggest a power-law decay of the spanning probability with L Z , where the exponent increases with the aspect ratio r. For fixed L X , it decays exponentially with r (see inset of Fig. 4).
It is possible to establish rigorously the off-critical power-law decay of Π(p), modifying the argument for the existence of anisotropic clusters presented above. Specifically, we can show that Π(p) ≥ L −θ X , where θ = θ(p, r) > 0, for any fixed r > 0 and p ∈ (p 2D , p c ). For that, let B be the diagonal band {(x, y); |x−L X /2| ≤ αn log(L X /n), |x − y| < 2n} in the center of the (x, y)-face of the cube, where α and n are constants setting the length and width, respectively. Let us say that the event B occurred if B is free of holes in the (x, y) plane i.e. if no sites in B are drilled in the z-direction. Also say that the event C occurred if there exists a path σ starting at height z = 0 and finishing at height z = L Z whose projection into the (x, y)-plane is contained in B and whose projection into the (x, z)-plane ((y, z)-plane) consists of sites that have not been drilled in the y-direction (x-direction). As discussed in detail in the Supplemental Material [7], for well chosen values of α and n the probability of the event C is bounded from below by a constant not depending on L X . Furthermore B and C are independent events. Since their joint occurrence implies the existence of a cluster including sites from the bottom and the top of the lattice, we conclude that where c 2 , c 3 , and θ are positive constants that depend on p.
The above argument also shows the existence of anisotropic clusters, sharpening the numerical results presented before. For p < p 2D , one has Π(p) ∼ e −c4L X , similarly to what happens for uncorrelated random percolation, where c 4 also depends on p. This is due to the fact that the projection of a path spanning the lattice into at least one of the coordinate planes is a path that spans the corresponding face, which has exponentially small probability in L X , due to the classical exponential decay of connectivity in the subcritical phase [26,27].
Conclusion. We find unexpected critical behavior when sequentially drilling holes through a solid cube until it is completely fragmented. At the critical density of drilled holes, a continuous transition is observed in a different universality class than the one of random percolation. We also numerically observe off-critical scale-free behavior that we can justify for a wide range of densities of holes using rigorous arguments. This model is a representative of more complex percolation models where sites are removed in a strongly correlated manner [28][29][30]. Examples are models where the set of removed sites is given by randomized trajectories, such as the so called Pacman and interlacement percolation models proposed to study the relaxation at the glass transition [31], enzyme gel degradation [32] and corrosion [33,34], as well as percolation models for distributed computation [35,36]. Other examples are percolation models explicitly introduce strong directional correlations as in the removal of cylinders [37] and different variants of the four-vertex model [38]. It would be interesting to explore up to which degree these models are in the same universality class or share common features.
While the fractal dimension of the largest fragment is consistent with the one of random percolation, all the other critical exponents are different. This has practical implications as the connectivity and transport properties do change considerably close to the threshold of connectivity. For example, we find the exponent of the order parameter to be substantially larger than for usual percolation which implies that the drilling transition is less abrupt. Since sites are removed along a line, it is necessary to remove more sites to produce the same effect in the largest fragment. We also find that, compared to usual percolation, the fractal dimension of the backbone is larger and the one of the shortest path is smaller, corresponding to a more compact backbone and therefore enhanced conductivity properties.
We acknowledge financial support from the ETH Risk Center, the To determine the position of the drilling transition, we consider in Fig. S1 the finite-size behavior of five different estimators: • χ ∞ refers to the maximum in the standard deviation of the largest cluster size. This is the position of the peaks in the curves in Fig. S6.
• M 2 is the position of the peak of the second moment, excluding the contribution of the largest cluster, as plotted in Fig. S7.
• Π is the average value of p where the first cluster spanning the lattice in z-direction appears [2].
• J is the point where the largest change in size of the largest cluster occurs, for every given sample, averaged over all samples [3].

II. DETAILED MODEL DESCRIPTION
Consider a simple cubic lattice of linear size L with free boundary conditions. To explore different fractions of drilled holes, we use the following setup. In the initial configuration, there are no holes, such that all N = L 3 unit cells (sites) are present (occupied).
We consider the following process. On each face of the cube i = x, y, z a site x 0,i is selected uniformly at random among the unselected sites, the corresponding hole is drilled, and all occupied sites in the direction perpendicular to the face and originating from x 0,i are removed (labeled as unoccupied). Therefore, in this step all occupied sites with coordinates x k,i = x 0,i + ke i are removed, where e i is the unit vector perpendicular to the cube face and k = 0, . . . , L − 1, such that at most L sites can be drilled in one step. This step is executed 3L 2 times until all holes have been drilled and consequently all sites in the system are gone. If we denote by (1 − p i ) the fraction of drilled holes in face i, we can characterize the state of the process with the control parameter p = p x = p y = p z . A system is considered to be connected if there exists a path of occupied nearest neighbors connecting two opposing faces of the cube, say in z direction. Let us denote by Π(p) the probability that the system is connected. Then, in the thermodynamic limit, Π(p) is a step function with Π(1) = 1 and Π(0) = 0 and with a transition at p = p c [15]. In Fig. S2 we see the connection probability as a function of the control parameter p for different linear system sizes L. It can be observed that the curves intersect around p c = 0.6339 ± 0.0005 (see also Fig. S1). Given its error bars, this value of p is higher than the site percolation threshold for the square lattice p c = 0.59274621 ± 00000013 [1,16,17] and lower than the cubic root of the site percolation threshold of the simple-cubic lattice: [p c (s.c.)] 1/3 ≈ 0.67795 [6].

III. OBSERVABLES AS FUNCTION OF p
The connection probability Π(p) is measured as the ratio of samples that are connected at parameter value p and the total number of samples. For example, if all holes have been drilled, we have p = 0 and Π(0) = 0, and likewise if no hole has been drilled, Π(1) = 1 (see Fig. S2). I. Critical exponents. The first column shows the symbol of the exponent. The second column shows the results for the drilling transition). For comparison, literature values for two-and three-dimensional classical percolation are given in columns three and four [1,5]. Here, β/ν is the exponent related to the order parameter, see Eq. (S1). γ/ν is the exponent related to the second moment [see Eq. (S7)] of the cluster size distribution and to the susceptibility [see Eq. (S3)]. 1/ν is the inverse of the correlation length critical exponent, as determined from the finite size scaling of threshold estimators (see Fig. S1) [1,6]. β is the order parameter critical exponent, determined from the behavior P∞ ∼ (p − pc) β [see Fig. 2(a)] [1]. γ is the susceptibility critical exponent, determined from the behavior M 2 ∼ |p − pc| −γ (see Fig. 2) [1,6]. dRS is the fractal dimension of the red sites in the backbone of the largest cluster [6][7][8][9]. dSP is the fractal dimension of the shortest path in the largest cluster [10]. dBB is the fractal dimension of the backbone of the largest cluster [11][12][13][14]. τ is the exponent of the cluster size distribution (see Fig. S11  Connection probability Π as a function of the control parameter p. The systems consist of N = L 3 sites where the largest considered linear system size is L = 1024 and the smallest one is L = 8. Results have been averaged over at least 10 7 samples for the smallest lattice size and 10 4 samples for the largest one. Π is measured as the fraction of connected samples. The curves for large lattices cross around pc = 0.6339 ± 0.0005 and Π(pc) = 0.21 ± 0.03. For clarity, the inset only shows the data for the four largest lattice sizes.
In the simple-cubic lattice, every site in the bulk has six nearest neighbor sites. If two occupied sites are nearest neighbors they are said to be connected. A set of connected occupied sites is called a cluster. Two distinct clusters are therefore separated by unoccupied sites. The number of sites forming the cluster is called its size. A single occupied site which has six unoccupied nearest neighbor sites is considered to be a cluster of unit size. The number of clusters per site C is defined as the ratio of the number of distinct clusters and the total number of sites N = L 3 . For critical percolation, this is known to be a lattice-dependent constant [21][22][23][24], the same seems to hold for the drilling transition, Fig. S3.
The largest cluster size per site P ∞ is defined as the ratio of the maximum cluster size among all clusters in the system and the total number of sites N = L 3 . For classical percolation, P ∞ acts as order parameter [1] which is zero in the disordered state and nonzero in the ordered one. The data for the drilling transition is shown in Fig. S4. For percolation close to the percolation thresh- old, one expects the following scaling behavior: where β is the critical exponent related to the order parameter, ν is the critical exponent related to the correlation length, p c is the percolation threshold, and F P is a scaling function [1]. If we compare classical percolation and drilling with respect to P ∞ against the fraction of removed cells, we observe that the cube gets disconnected shows the data for the four largest lattice sizes rescaled with pc = 0.6339, 1/ν = 0.915, and γ/ν = 2.0452. The ratio γ/ν is the one of three-dimensional percolation [6]. 1/ν = 0.915 ± 0.010 is consistent with the scaling behavior of P∞ (see Fig. S4) and the estimators for pc (see Fig. S1), but it differs from the value for percolation, 1/ν = 1.1450 ± 0.0007 [6,[18][19][20].
at lower fractions of removed cells in the case of classical percolation, see Fig. S5.
The standard deviation of the largest cluster size per site χ ∞ is defined as the square root of the variance of P ∞ : where · S indicates averaging over independent realizations (see Fig. S6). For percolation close to the percolation threshold, one expects the following scaling behavior: where d is the spatial dimension of the lattice, γ is the critical exponent related to the susceptibility, and F χ is a scaling function. The second moment of the cluster size distribution is defined as where the sum runs over all clusters in the system and s k is the size (number of sites) of cluster k. Excluding the contribution of the largest cluster size, we arrive at the following definition: where s max is the largest cluster size. In Fig. S7, we see Close to the percolation threshold, one expects the following scaling behavior: where F M2 is a scaling function. According to the scaling forms in Eq. (S1) and (S7), one expects to observe and at the threshold p = p c . The corresponding measurements are shown in Fig. S8.

IV. ALGORITHMS AND BOUNDARY CONDITIONS
To simulate the drilling model with an algorithm of complexity linear in the number of sites N = L 3 , we adopt the following strategy. We consider initially a system with all sites occupied. Starting from this setup, the holes are drilled according to the processes described above. For each site, the fraction of drilled holes (1 − p x ) at which it is removed is recorded. Once all holes have been drilled, we proceed in the opposite direction, and analyze the following percolation problem. Suppose we start from a lattice where all sites are unoccupied (p x = 0). Now, p x is incremented and the sites of the lattice become occupied in the inverse order of the one in which they have been drilled. This procedure allows to keep track of the properties of the clusters of occupied sites as function of the control parameter p x [2,16,17]. Random numbers have been generated with the algorithms discussed in Refs. [29][30][31]. As a further test, we also considered alternative boundary conditions, obtaining results consistent with the ones obtained for free boundary conditions. In particular, we employed a single-cluster growth method, similar to the Leath algorithm for classical percolation [32,33]. This method works at a fixed occupation probability p. We begin with a lattice of length L equal to a small multiple of three. The holes in this initial cube are present with probability 1 − p and absent with probability p. If the site in the center of the lattice is occupied (this is the case with probability p 3 ), a cluster is grown starting from there, connecting to all occupied neighbors. In case the cluster touches a site on the boundary of the lattice, the linear lattice size is increased by a factor of three, the information on the drilled holes is accordingly propagated to the new lattice and the cluster growth continues. At a certain cutoff lattice size L, the cluster growth is terminated. This measurement gives the number of clusters as function of the size of the cluster at the origin, as shown in Fig. S10, at p = p c . Because the probability of a randomly picked site to be in a cluster of size s is proportional to s, this measures a size distribution which behaves as sp(s). The corresponding power-law exponent is measured as τ − 1 = 1.165 ± 0.01, see Fig. S11. The influence of varying p on the cluster size distribution is shown in Fig. S12. Finally we verified that the cutoff of the cluster size distribution, s cut , scales as expected from the fractal dimension of the largest cluster: We measured the fraction f (L) of clusters grown from the origin that reach the boundary of a box of length L and found that f (L) ∼ L −0.48±0.10 (not shown). One expects f to be Probability ∏(p=0.6339) S9. (a) shows the connection probability Π as function of the lattice length LZ , with p = 0.6339, for different lattice aspect ratios r. Π saturates to a constant for large LZ . As can be seen in (b), for large aspect ratios r, the connection probability decays exponentially in r, similarly to classical percolation [25][26][27][28]. The lattice size is LX × LY × LZ with LX = LY and LZ = rLX and the connection probability is measured in z-direction. Results are based on at least 10 8 samples.
related to the cluster size distribution p(s) ∼ s −τ by assuming s cut ∼ L d f [34,35] and the validity of the scaling relation τ = 1 + d/d f [1,36], in analogy to classical percolation. In Fig. S8, we measure d f −d = −0.50±0.04, in agreement with the behavior of f .

V. CLUSTER SHAPE AND POWER-LAW DECAY
To study the geometry of the drilling transition clusters, we considered the eigenvalues and eigenvectors of the cluster inertia tensors [37][38][39], see Fig. S13 to S15. The results show that the clusters of the drilling transition are more anisotropic and aligned with the cube edges as the ones of classical 3D percolation.
We now present the details involved in the proof of the polynomial decay for the spanning probability Π(p) as a function of the lattice length L X (recall that L Y = L X and L Z = rL X with a fixed aspect ratio r). We fix a p ∈ (p 2D , p c ). We consider the diagonal band contained in the (x, y)-plane B = {(x, y); |x−L X /2| ≤ αn log(L X /n), |x− y| ≤ 2n} and denote byB = {(x, y, z); (x, y) ∈ B, 0 ≤ z ≤ L Z }. Recall the definition of the event C that there exists a path σ contained inB starting at height z = 0 and finishing at height z = L Z consisting of sites whose projections onto the (x, y) and (y, z)-planes are free of holes. Our goal is to show that if α and n are chosen sufficiently large then C has probability strictly larger than 0 uniformly in L X , i.e. that there exists a δ > 0 depending on p and r only, such that P [C] ≥ δ. The argument involves a one-step renormalization. Let us consider the rectangles B X = [0, αn log(L x /n)] × [0, L Z ] and B Y = [0, αn log(L Y /n)] × [0, L Z ] contained in the (x, z) and (y, z)-planes respectively. We tileB with cubes of side length n, inducing a tiling of B X and B Y with squares of side length n L X . This allows us to see B X (respectively B Y ) as a rectangle in a renomalized square lattice whose sites correspond to n × n squares composed of n 2 sites of the original coordinate (x, z) and (y, z) planes respectively. For a given n × n square, we say that it is occupied if it is surrounded by a circuit of sites free of holes contained in its eight neighboring n × n squares (where two squares are considered neighbors if they intersect each other). Since p > p 2D the probability that a square is occupied can be made arbitrarily high by choosing n large enough. Thus we can choose n and α large enough so that the probability of finding a path of cubes inside B traversing it from bottom to top and whose projections are occupied squares in B X and B Y is strictly larger than 0 uniformly in L X . Note that the fact that the projections are occupied squares, and the fact that the circuits of neighboring occupied squares intersect each other, assures that one can find a path inside B like in the definition of the event C. (a) Histogram of the ratio of the smallest over the largest eigenvalue of the inertia matrix [40] for the drilling model (p = pc = 0.6339) and classical percolation (p = pc = 0.3116077 [6]) clusters in their center of mass frame, measured for different lattice sizes L, for all clusters with size at least 10. The inertia matrix consists of the following elements: I ik = (x 2 l δ ik − xix k ), where the sums are over all sites in the cluster. The different behavior of the left shoulders of the histograms indicates the relatively large amount of clusters with elongated shape in the drilling model.   S14. (a) Histogram of the maximum alignment, with a coordinate axis, of the eigenvector of the inertia matrix with the smallest eigenvalue at criticality. We determine the eigenvector of the inertia matrix corresponding to the smallest eigenvalue. For each of the three coordinate axes, the absolute value of the cosine of the angle between the vector and the axis, ai, is calculated. Then, the maximum alignment is defined as the highest ai for i = x, y, z. For the drilling model, the distribution of the maximum alignments shows a peak close to unity, corresponding to clusters whose anisotropy is aligned with the coordinate axes. (b) Histogram of alignments for the drilling model, for different values of the control parameter p (and lattice size L = 256).  Fig. S14, however considering the alignment of the eigenvector with the lowest eigenvalue with a fixed coordinate axis (say, the x-axis, ax). While the distribution of the alignment seems to become nearly uniform for classical percolation, for the drilling model one observes peaks close to zero and unity, compatible with the anisotropy of the clusters being aligned with the coordinate axes. (b) Alignment histograms as in (a), for the drilling model at different values of p (and with lattice size L = 256).