Online SI : Collective ordering of colloids in grafted polymer layers

with rb the blob radius. In this regime the probability of 3 blobs overlapping is sufficiently small, we neglect many body effects and interactions between blobs are pairwise additive. By making blobs smaller we need more of them to describe a single polymer and consequently computational cost increases. Blob size was chosen as rb = σ/6 with σ the colloid diameter, which was proven to be a good tradeoff between accuracy and efficiency. With this chosen blob size the effective 3D blob density in our simulations always stays within limitations of the model (Eq. 1), i.e. the effective blob density on snapshots shown in Figure 6 of the main text is about ρblob ≈ 1/8r b . In order to additionally verify that the blob picture is appropriate, we performed insertion free energy calculations for two different blob sizes rb = σ/6 and r ′ b = σ/12, the results were essentially the same. We also performed many colloid simulations for ρ = 3, g′ = 6, with different ratios of blob sizes to colloid diameter of (2, 3, 4), again, the resulting morphologies and their typical width were essentially the same. With this model we performed canonical and grand canonical Monte Carlo simulations. Surface area, number of chains and number of blobs in a chain remained constant throughout a single run. Brush conformation is changed via a single blob translational moves. Colloids are also moved around with single translational moves or (in a Grand-canonical simulation) randomly inserted or deleted. Simulation box size was always 20 times 20 colloid diameters in x and y with periodic boundary conditions. Box size in z direction was much grater than the brush thickness so that none of the particles ever interacted with the top surface of the box. Single colloid insertion free energy Fp(z) was calculated using Wang-Landau technique [2, 3] biasing in the z coordinate of the colloid with the biasing function ψ(z). The bin size was chosen to be σ/12. Initial value of the modification factor was f = 1, when the ratio between the maximum and minimum value of the histogram of all visited states became maxH(z) minH(z) < 1.5 the new value for the modification factor was assigned f = f/4. This scheme was continued until the modification factor became sufficiently small f < 10−5, after which the WangLandau algorithm was terminated (setting f = 0) and a long simulation using the obtained biasing function was performed. The final free energy was calculated as βFp(z) = ψ(z)− ln (H(z)) + C, with H(z) the final histogram of visited states, β = 1/kBT and C a constant. C was determined such that the free energy is zero when the colloid is far above the brush and not interacting with it Fp(z0) ≡ 0, with z0 hb, hb being the brush height. For a given brush, with certain ρ and lp, the free energy calculations were done and averaged over N = 50 parallel systems, each system having a different random realisation of the chain anchoring points

Polymers are coarse-grained using a chain of soft repulsive blobs.Ideally the size of a single blob should be much smaller than the colloid diameter so that from the blob's point of view the colloid looks almost like a flat wall.We additionally need to take into account that the blob model used is valid for densities smaller than approximately one blob per blob volume [1] with r b the blob radius.In this regime the probability of 3 blobs overlapping is sufficiently small, we neglect many body effects and interactions between blobs are pairwise additive.By making blobs smaller we need more of them to describe a single polymer and consequently computational cost increases.Blob size was chosen as r b = σ/6 with σ the colloid diameter, which was proven to be a good tradeoff between accuracy and efficiency.With this chosen blob size the effective 3D blob density in our simulations always stays within limitations of the model (Eq.1), i.e. the effective blob density on snapshots shown in Figure 6 of the main text is about ρ blob ≈ 1/8r 3 b .In order to additionally verify that the blob picture is appropriate, we performed insertion free energy calculations for two different blob sizes r b = σ/6 and r b = σ/12, the results were essentially the same.We also performed many colloid simulations for ρ = 3, g = 6, with different ratios of blob sizes to colloid diameter of (2,3,4), again, the resulting morphologies and their typical width were essentially the same.
With this model we performed canonical and grand canonical Monte Carlo simulations.Surface area, number of chains and number of blobs in a chain remained constant throughout a single run.Brush conformation is changed via a single blob translational moves.Colloids are also moved around with single translational moves or (in a Grand-canonical simulation) randomly inserted or deleted.Simulation box size was always 20 times 20 colloid diameters in x and y with periodic boundary conditions.Box size in z direction was much grater than the brush thickness so that none of the particles ever interacted with the top surface of the box.
Single colloid insertion free energy F p (z) was calculated using Wang-Landau technique [2,3] biasing in the z coordinate of the colloid with the biasing function ψ(z).The bin size was chosen to be σ/12.Initial value of the modification factor was f = 1, when the ratio between the maximum and minimum value of the histogram of all visited states became max H(z) min H(z) < 1.5 the new value for the modification factor was assigned f n = f o /4.This scheme was continued until the modification factor became sufficiently small f < 10 −5 , after which the Wang-Landau algorithm was terminated (setting f = 0) and a long simulation using the obtained biasing function was performed.The final free energy was calculated as βF p (z) = ψ(z) − ln (H(z)) + C, with H(z) the final histogram of visited states, β = 1/k B T and C a constant.
C was determined such that the free energy is zero when the colloid is far above the brush and not interacting with it F p (z 0 ) ≡ 0, with z 0 h b , h b being the brush height.For a given brush, with certain ρ and l p , the free energy calculations were done and averaged over N = 50 parallel systems, each system having a different random realisation of the chain anchoring points Total insertion free energy of a single colloid in an external field is obtained by summing up the contributions from polymer repulsion and external field, since polymer chains are not affected by external field and the two terms are decoupled F (z, g ) = F p (z) + g z/(βσ). (3)

II. 2D HOLE-SIZE DISTRIBUTION
Here we derive the analytical expressions for the distribution of the two-dimensional circular voids on surfaces characterized by quenched disordered (random Poissonian) or ordered (square crystal) distribution of the anchoring points.

Random anchoring points.
Fixing the frame of reference to an arbitrary anchoring point, we want to evaluate the probability p + (r) to find another point at a distance between r and r + dr from it (with no points between 0 and r).On a randomly grafted surface there are no correlations and we have where ρ is the mean grafting density.'+' designates probability 'to find' and '−' probability 'not to find'.The probability not to find an anchoring point in a circle with radius r e −p + (r ) = e − 2πr ρdr = e −πr 2 ρ .
(5) So, in total, the probability to find a first anchoring point a distance r from a randomly chosen point is We would get the same result by differentiating p +1 (r) = −∂(P − (0, r))/∂r.
The spacing a between the points on a square lattice is constant, given by the mean density, ρ = 1/a 2 .The probability to find an anchoring point a distance r from a randomly chosen point is 2πrdr for r ≤ a/2.Then p +1 (r) = A 2πrdr, where A is a normalization factor.When a/2 < r ≤ a √ 2/2, we have to work out the arc length l of a circle with radius r that falls inside a square with the side a.A straightforward derivation gives l = 2r(π − 4 arccos a 2r ) and so p +1 (r) = A 2r π − 4 arccos a 2r dr.Using ρ = 1/a 2 and normalizing we get A = ρ, therefore in total (7) These analytically derived expressions are displayed in the Figure 1b) in the main text of the manuscript together with the numerically obtained results.

III. INSERTION FREE ENERGY SCALING
Insertion free energies for disordered grafting at various chain lengths and grafting densities are shown on Figure 1 .We immediately notice two general scaling relations.The brush height h b scales approximately linearly with chain length l p and the insertion free energy increases approximately linearly with the grafting density (also see Figure 2 in the main text).We also observe that the brush height increases slightly with increasing ρ.
In the "mushroom" regime the height of the polymer layer depends on R g of individual chains which, for selfavoiding walk polymers, scales as R g ∝ l 0.58 p , it does not depend on ρ.On the other hand in the "brush" regime monomer crowding stretches the polymers and brush height also depends on ρ.For self-avoiding chains with the exponent ν = 0.35 [4].Since we observe a dependence of the brush height on ρ we are not in the "mushroom" regime.Therefore we can assume that we are in the "brush" regime and approximately write that changing the grafting density from ρ to ρ * , while keeping the chain length l p and colloid diameter σ constant, changes the insertion free energy where we have to rescale the height of the free energy barrier and also include the brush height dependance on ρ.
We use a general notation that the insertion free energy depends on the colloid diameter, grafting density, chain length and vertical position F = F (σ, ρ, l p , z).Using a different sized colloid σ * effectively means we change our unit of length For long chains, where brush height is larger then about 4 colloid diameters, we also observe (on Figure 1) that increasing the chain length does not change the barrier height, it merely stretches the free energy profile, so we can approximate that Collecting all scaling relations together (Eq.9, 10, 11) we get a master scaling equation which should work for high enough grafting densities and chain lengths.Lower right corner of Figure 1 shows the scaling performance.Ideally all curves should fall on a single master curve, we see that the scaling works rather well when the grafting density and chain length are large enough.The blue curve that slightly sticks out is the one at ρ = 1 and l p = 40 (corresponding to ρ ≈ 2), where we are no longer in the "brush" regime, but somewhere between "brush" and "mushroom".
The barrier height scales linearly with the grafting density (see Fig. 2e) of the main text) -the ideal gas like behavior -here we attempt to provide a hand waiving explanation for that behavior.For the range of grafting densities studied we are above the "mushroom" regime and there is appreciable contribution to the pressure (and insertion free energy) from the second virial coefficient (due to chain -chain repulsion), so pressure does not scale linearly with 3D blob density but with an exponent larger than 1.On the other hand at higher grafting densities (brush regime) the polymer layer height increases with the grafting density, therefore ρ 3D b does not scale linearly with the grafting density but with an exponent lower than 1.The two effects approximately cancel out and we find that the insertion free energy scales approximately linearly with the grafting density.

IV. DETERMINING κ FROM SIMULATIONS
As explained in the main text the balance of hydrostatic and polymer pressure gives with κ a constant.The value of the critical density should therefore linearly depend on the brush height and strength of external force field, conversely it should not depend on the grafting density.Parameter κ is determined from the plot of ρ 0 as a function of g h on Figure 2.For a given l p and g a set of systems at various ρ was considered, then ρ 0 was determined from fitting a linear function to plots of η vs. ρ (like the one shown in the inset of Figure 2).For each system, η was determined from a 2D blob density plot (with mesh size a m = σ/3) by determining the fraction of low density area with ρ m < ρ cut m , ρ m being the local density in a particular mesh point.The cutoff density was chosen as ρ cut m = ρ high m /5, ρ high m being the average blob density in the 'high' density area.That was found to give a good distinction between low and high density areas, but the results are not sensitive to a particular choice ρ cut m as long as it is reasonable.We then also have a relationship ρ 0 ρ high m /l p .Brush monomer density is a decaying function of z (see Figure 2 of the main text) so there are many ways in which to determine the brush height.Our choice was to determine h b by considering the blob height distribution For each set of systems at given lp and g and various ρ we determine h b with Eq. 14 and ρ0 with Eq. 13.Horizontal lines depict errors.A Linear function is fitted to the data and its slope corresponds to κ = 0.11.Our theory says that the brush height should be independent of ρ because colloids always compress the brush until the effective brush density reaches ρ0.We see this assumption is quite accurate, it only deviates for η very close to 1 (data marked with green circles), where there is considerable amount of polymer mass being squashed beneath the colloids.Inset: Dependence of η on ρ for random (red circles) and square (black squares) anchoring with lp = 40 and g = 6, this data corresponds to snapshots on Fig. 6 of the main text.
meaning that polymer layer height is defined where there is 95% of polymer mass under that height.This definition gives us a layer height that is roughly equal to the average height of a single colloid sitting on the top of the layer (e.g.see ρ = 4 curve on Figure 3 of the main text).

V. DETERMINING Cex AND Csq
Here we attempt to determine the value of the constant C ex which specifies the average free energy penalty per unit volume per grafting density for inserting an object into a polymer layer C ex = F ex /((V + V γ )ρ), where V is volume of the object and V γ is volume of the depletion layer.We can approximately determine C ex by using the single colloid insertion free energy data (Figure 1).We can compute the average insertion free energy over the polymer layer by integrating where h b is the brush height determined by Eq. 14.We plot F (ρ) for different chain lengths on Figure 3. Fitting a linear function to the data we find that the average insertion free energy (per colloid) per grafting density C * ex ≈ 6k B T σ 2 .Then we can compute Using γ ≈ 0.2 and (V col + V γ) = 4π(r + γ) 3 /3 = 1.43 .
The constant C sq specifies the Free Energy penalty for squashing a unit length of a chain.We assign a 1k B T penalty for every deGennes blob (with diameter ξ) being squashed, so C sq = k B T /ξ.Using ξ = 1/ √ ρ 0 we get For data and snapshots presented on Figure 6 of the main text ρ 0 = 6 and C sq ≈ 2.4k B T

FIG. 2 .
FIG.2.Determining κ from simulation snapshots.For each set of systems at given lp and g and various ρ we determine h b with Eq. 14 and ρ0 with Eq. 13.Horizontal lines depict errors.A Linear function is fitted to the data and its slope corresponds to κ = 0.11.Our theory says that the brush height should be independent of ρ because colloids always compress the brush until the effective brush density reaches ρ0.We see this assumption is quite accurate, it only deviates for η very close to 1 (data marked with green circles), where there is considerable amount of polymer mass being squashed beneath the colloids.Inset: Dependence of η on ρ for random (red circles) and square (black squares) anchoring with lp = 40 and g = 6, this data corresponds to snapshots on Fig.6of the main text.
p b (z), we define h b by an integral h b 0 p b (z)dz ≡ 0.95,

ElectronicFIG. 3 .
FIG. 3. Average insertion free energies for disordered grafting.A linear function is fitted to the data at lp = 40 and its slope corresponds to C * ex = 6kBT σ 2 .