Debye vs . Casimir: controlling the structure of charged nanoparticles deposited on a substrate †

Fine-tuning the interactions between particles can allow one to steer their collective behaviour and structure. A convenient way to achieve this is to use solvent criticality to control attraction, via critical Casimir forces, and to control repulsion via the Debye screening of electrostatic interactions. Herein, we develop a multiscale simulation framework and a method for controlled deposition of quantum dots to investigate how these interactions a ﬀ ect the structure of charged nanoparticles deposited on a substrate, altogether immersed in a binary liquid mixture intermixed with salt. We consider nanoparticles and substrates favouring the same component of the mixture and ﬁ nd that the critical Casimir interactions between the nano-particles become drastically reduced at the substrate. In particular, the interactions can become a few k B T weaker and their decay length a few orders of magnitude smaller than in the bulk. At o ﬀ -critical compositions, the decay length increases upon approaching criticality, as expected, but the interaction strength decreases. With molecular dynamics simulations and experiments, we reveal that the nanoparticles can self-assemble into crystalline clusters which form superstructures resembling cluster ﬂ uids and spinodal morphology. The simulations additionally predict the formation of fractal-like nanoparticle gels and bicontinuous phases. Our results demonstrate that charged nanoparticles in a salty binary liquid mixture provide exciting opportunities to study the formation of complex structures experimentally and theoretically, which may lead to applications in optoelectronics and photonics.


Introduction
The ability to tune the interactions between particles provides a means to control their collective behaviour and structure.Perhaps the simplest way to tune the range and the strength of interactions between charged particles is to vary the concentration of ions in solution.In the case of identical particles, the interactions are repulsive and decay exponentially with the Debye screening length, which is proportional to inverse of the square root of the ion concentration.The Debye length can be well-controlled in the range from a nanometre to a hundred of nanometres, 1 corresponding to the ion concentration ranging from about 1 M to 1 μM, respectively, for water-2,6-lutidine mixture at room temperature.
Another way to control interactions among particles is through immersion in a near-critical fluid.][4][5][6][7][8][9][10][11][12][13] These forces are called critical Casimir forces, in analogy with the Casimir forces due to quantum fluctuations, 14 as they arise due to fluctuations of the order parameter associated with the solvent criticality; 15 for a binary liquid mixture, for instance, such an order parameter is the difference in the concentrations of the two components.Critical Casimir forces are attractive for identical particles, and their range and strength can be tuned by changing the composition of the fluid or by varying temperature. 137][18] Such micrometer-sized colloids can be observed and studied in realtime using confocal laser scanning optical microscopy, which is inapplicable to colloids smaller than the diffraction limit.However, fine temperature control of the order of few milli Kelvins is needed to reach correlation lengths ξ comparable to the particle size, in order for the critical Casimir interactions among the colloids to become significantly strong.The maximum achievable ξ reported so far is about 200 nm (see, e.g., ref. 19) so that the colloid-colloid interactions extend to about one colloidal radius at most.
Using nano-sized particles effectively extends the range of critical Casimir interactions up to a few particle radii, while simultaneously requiring a coarser tuning of the temperature to control the interaction strength and range.1][22][23][24][25] For instance, Hopkins et al. used classical density functional theory to calculate the effective potential between nanoparticles dispersed in a binary liquid mixture, 20 and El Fassi et al. calculated the critical Casimir force between nanoparticles immersed in a crosslinked polymer blend. 21Lu et al. studied glass transitions in a suspension of silica nanoparticles in a water-lutidine mixture near its consolute point, 22 while Marino et al. used the critical Casimir forces to assemble and control quantum dot superstructures. 23,24erein, we study how the combination of electrostatic repulsion, via Debye screening, and critical Casimir attraction determines the structure of charged nanoparticles depositing on a substrate.We consider substrates and nanoparticles with the same preference for the fluid phase or component of a binary liquid mixture (e.g., hydrophilic nanoparticles on hydrophilic substrates) and investigate how nanoparticles structure depends on various parameters, such as the composition of a mixture, temperature and salt concentration.To this end, we develop a multiscale simulation framework.We first employ Monte Carlo simulations of the Ising model, which mimics a binary liquid mixture or a simple fluid close to criticality, to obtain an effective critical Casimir potential between the deposited nanoparticles.We then use this potential in molecular dynamics simulations to study their collective behaviour.To compare the simulated structures to experiment, we develop a novel method to drive the deposition of quantum dots on silicon substrates that does not rely on DNA functionalisation 26 or nanoscale patterning. 27The deposited superstructures reproduce the main features of the simulations well, demonstrating a novel pathway to quantum dot epitaxy through critical Casimir forces.

Interaction of nanoparticles at a substrate
Binary liquid mixture, Ising model and thermodynamic paths Our goal in this article is to understand the behaviour of nanoparticles at a substrate, as shown in Fig. 1a, driven by the interplay between the Debye-screened electrostatic repulsion and critical Casimir attraction.The latter emerges in solvents such as binary liquid mixtures in the vicinity of their critical points.A computationally feasible and reliable method to study such interactions is to employ the Ising model (eqn (6) in Methods), which belongs to the same universality class as a binary liquid mixture. 28Magnetization m in the Ising model corresponds to concentration c of a binary mixture and bulk field h to the difference of the chemical potentials Δμ of the two species of the mixture (for the details of this correspondence, see ref. 29-31).Schematic phase diagrams of the two systems are drawn in Fig. 1b and c.
In experiments with binary liquid mixtures, one typically keeps the concentration fixed and varies the temperature, as denoted schematically by a thin vertical line in Fig. 1b.In simulations of the Ising model, however, it is more convenient to vary either the bulk field or temperature, while keeping one of these variables constant (the lower plot in Fig. 1c).This has been done in virtually all theoretical and simulation studies on critical Casimir forces [10][11][12][13][32][33][34] (see, however, ref. 35 and 36 where such forces were computed in the canonical ensemble). In ths work, in order to mimic our experimental system as closely as possible, we chose to follow the path of constant bulk magnetization, which corresponds to constant concentration of a binary liquid mixture.To this end, we first performed Monte Carlo (MC) simulations of the Ising model in bulk and determined the bulk magnetization m for various values of the bulk field h and reduced temperature τ = (T − T c )/ T c , where T c is the critical temperature.We then selected the (h, τ) pairs that provide the same bulk magnetization (Fig. S1 and Table S1 †) and performed MC simulations for two nanoparticles at a substrate and in bulk for these values of h and τ.In addition to these simulations, for comparison we also carried out simulations along the conventional paths of constant h and τ.

Critical Casimir interactions
Our MC simulations reveal that the presence of substrate can significantly modify the critical Casimir interactions between nanoparticles.Fig. 1d compares the interaction potentials between two identical particles in bulk and at a substrate, where the substrate favours the same fluid component as the particles, for the parameters close to the critical point.At the substrate, both the depth and the extent of interactions are significantly reduced.This is because the substrate promotes an environment (i.e., attracts a fluid component) preferred by the particles (Fig. 1e), which reduces the deformations of the magnetization (or concentration) profile between the particles and consequently the strength of the interaction.
The effect of substrate becomes less pronounced as one moves away from the critical point, that is, upon increasing h and τ (Fig. 1f ).In the considered case, the bulk field in the Ising model favours magnetization opposite in sign to what is preferred by the nanoparticles.This leads to the formation of a capillary-like bridge connecting the two particles (Fig. 1g) and correspondingly strong but short-range forces. 34,37,38Since the bridge formation is similar at the substrate and in bulk, the interaction potential changes little when the particles are brought to the substrate.
To quantify the reduction of the critical Casimir interactions by substrate, we extracted the decay lengths λ C of the critical Casimir potentials by fitting the MC data to an exponentially decaying function of separation where β = (k B T ) −1 , with k B the Boltzmann constant and T temperature; D = r − 2a is the surface-to-surface distance between the particles, and b and A are two additional fitting parameters.Here the 1/D dependence (for D ≫ b) is related to the behaviour of the potential in the vicinity of criticality, 4,5,39 b describes the non-universal behaviour at short separations and A is the amplitude.This function fits the MC data in the whole range of distances remarkably well (solid lines in Fig. 1d and f ).Fig. 1h shows that the substrate reduces the decay length λ C by a few orders of magnitude in the region where λ C in bulk is extremely large (>10a), but it remains practically unaltered in the region where λ C ≲ a.The Casimir decay length λ C (for particles at the substrate) is shown as a heat-map in Fig. 1i along three different paths.Although reduced as compared to the bulk value, λ C is still large close to the critical point (τ = 0, h = 0) but decreases as τ or h increases.For constant bulk field, the decay length varies little with τ and decreases significantly only sufficiently far from τ = 0.

Interplay of the Debye-screened repulsion and critical Casimir attraction
The total interaction potential between the deposited, charged nanoparticles consists of the above-discussed critical Casimir attraction and electrostatic and Lennard-Jones interactions.However, we consider nanoparticles with the dielectric constant close to the solvent's dielectric constant, implying small or vanishing van der Waals (dispersion) forces.This assumption is consistent with our experiments, where the relative dielectric constant of nanoparticles (quantum dots) ε ≈ 5-6.25, 40,41 while ε ≈ 7 for the critical water-lutidine mixture. 42For non-vanishing dispersion forces, the electrostatic repulsion (see below) will have to be re-adjusted to compensate partially for the der Waals attraction.We note that both van der Waals and electrostatic interactions change very little with temperature in the range of temperatures considered (20 °C to 33 °C), 43,44 compared to drastic changes of the critical Casimir interactions. 45Thus, the system behaviour can be steered by the interplay between the Debye-screened repulsion, through varying salt concentration, and critical Casimir attraction, via fine temperature control.
In this work, therefore, we neglect dispersion interactions and combine the computed critical Casimir potentials with the Debye-screened electrostatic potential, taken in the DLVO form, and short-ranged Lennard-Jones repulsive potential (eqn ( 4) and (3) in Methods).Fig. 2a demonstrates how this combination creates a minimum (−E depth ) and an energy barrier (E max ) in the total interaction potential.We calculated E depth and E max along the thermodynamic paths used in our MC simulations (Fig. 1i) at different values of the Debye screening length λ D as well as for varying λ D at fixed values of bulk field h and reduced temperature τ.This is shown in Fig. 2b and c in the plane of λ D and the extracted decay length λ C of the critical Casimir potential at the substrate.
Perhaps surprisingly, the depth of the potential decreases with increasing λ C , that is, the interactions become weaker upon approaching criticality (Fig. 2b).This is in contrast to the behaviour at the critical composition in bulk, where both the potential depth and the extent of interactions increase upon approaching the critical temperature. 6,13The substrate screens the interactions in the immediate vicinity of criticality, reducing the potential depth significantly.The increase of E depth with moving away from the criticality can be attributed to the formation of capillary bridges 34,37,38 (Fig. 1g).We note that we observe a similar behaviour at off-critical compositions also in bulk, but its magnitude is less pronounced (Fig. S2 †).
A large Debye screening length λ D can lead to the appearance of an energy barrier E max and hence to the repulsive tail in the total interaction potential (Fig. 2a and c).7][48][49][50] We will see that the emergence of such repulsive tails in our potentials can influence the nanoparticle behaviour in a nontrivial way.

Collective behaviour of nanoparticles at the substrate
To investigate the collective behaviour of nanoparticles at the substrate, we performed 2D molecular dynamics (MD) simulations using the total interaction potentials discussed above.To characterize the obtained structures, we calculated radial distribution functions (RDFs), which describe how density varies as a function of distance from a reference particle.In addition, we computed distributions of clusters formed by the nanoparticles and the average cluster sizes.A nanoparticle was considered part of a cluster when the distance to at least one of the cluster's nanoparticles was smaller than 2.5 the particle radius.To gain insights into the structure within a cluster, we calculated the hexatic order parameter ψ 6 for the nanoparticles belonging to this cluster, where ψ 6 = 1 for the perfect hexagonal lattice (eqn (7) in Methods).

Dependence on the decay length of the critical Casimir interactions
To study the dependence of the nanoparticle structure on the decay length λ C of the critical Casimir interactions, we chose the thermodynamic path of constant (non-zero) magnetization in the Ising model, which corresponds to a constant (off-critical) concentration in a binary liquid mixture (Fig. 1).To single out the effect of λ C , we took a low value of the Debye screening length, λ D /a = 1, so that the long distance behaviour is mainly determined by the critical Casimir interactions.
Fig. 3a shows that as the decay length λ C decreases, the depth of the total interaction potential increases (see also Fig. 2a and b).For a small λ C ≈ 0.53a, the potential depth reaches the value of about 8k B T. Such strong attraction leads to the formation of clusters of nanoparticles with almost perfect hexagonal ordering within the clusters.This is demonstrated by the short-distance behaviour of RDF and by the order parameter ψ 6 ≈ 1 (Fig. 3b and d).Remarkably, the nanoparticle clusters form a structure redolent of colloidal gels, consisting of a network of fractal-like aggregates. 51,52s the decay length λ C increases, the clusters coarsen and the structure inside the clusters becomes more amorphous.This is manifested by the decreasing order parameter ψ 6 and occurs because the depth of the potential decreases with increasing λ C .The nanoparticles can thus detach from the clusters more easily, which leads to a decreasing average cluster size.For the largest decay length studied (λ C ≈ 3a), the nanoparticles assemble into a structure reminiscent of the spinodal morphology.Similar coarsening of nanoparticle clusters with approaching criticality can be obtained at the critical composition (Fig. S3 †) and along the thermodynamic paths of constant reduced temperature τ or constant bulk field h, pro-

Dependence on the Debye screening length
The strength and the range of electrostatic repulsion can be controlled by the Debye screening length λ D through varying the salt concentration (eqn (4) in Methods).To study how the nanoparticle structure is affected by λ D , we chose the parameters providing a large value of the decay length λ C of the critical Casimir attraction, so that we could obtain a transformation between the potentials with attractive and repulsive tails.
Fig. 4 shows that the total interaction potential and the structure change significantly upon varying λ D .At a high salt concentration, implying a small λ D , the repulsive electrostatic interactions are strongly screened, and the interaction between the nanoparticles is attractive at all separations (λ D /a = 0.5 in Fig. 4a).This leads to the formation of coarsened structures, consisting of randomly arranged, partially connected clusters, within which the nanoparticles are hexagonally ordered.Again, the order is manifested by well-developed peaks in the RDF and by a high value of the hexatic order parameter (Fig. 4b and d).With increasing λ D , i.e., with reducing the salt concentration, a small energy barrier develops in the interaction potential, which impedes the nanoparticles from coming close to each other.This breaks the hexagonal order and we observe the growth of spinodal-like structures, characterized by smaller cluster sizes and lower values of the hexatic order parameter (Fig. 4d and e).Upon further increase of λ D , the nanoparticle clusters break and form a dispersed liquidlike structure without apparent long-ranged order.In this case, the RDF quickly decays to unity, and the hexatic order parameter and the average number of clusters decrease.The latter saturates at s ≈ 8 already for λ D /a ≳ 3. Interestingly, the hexatic parameter saturates at ψ 6 ≈ 0.4, which suggests that the particles tend to form local hexagonal order even in the dispersed state.

Effect of density
We have carried out MD simulations for various densities of the deposited nanoparticles for the parameters from Fig. 3 and 4. In most cases, we found that the density has a minor effect on the nanoparticles structure (Fig. S7-S10 †).For increasing density, the average cluster size increases, as small clusters unite to form bigger structures, but the ordering inside the clusters does not change dramatically.Perhaps not very surprisingly, our simulations show that the most interesting behaviour is obtained in the case of mermaid potentials with sufficiently large potential depths and relatively low energy barriers.One such case is shown in Fig. 5 and discussed below.
At a low surface area fraction (η = 0.1), the nanoparticles self-assemble into small polydisperse clusters that arrange themselves into a liquid-like structure with the tendency for ordering (Fig. 5b).This is manifested by long-ranged oscillations apparent in the RDF (the inset in Fig. 5c).However, we have not obtained any well-defined crystal structure of these clusters, likely due to their polydispersity (Fig. 5d).It is worth noting that similar cluster-fluid structures have been observed in experiments with quantum dots at a water-air interface. 53he cluster size distribution exhibits a peak that widens and moves toward larger s as η increases (Fig. 5d).Correspondingly, the average cluster size increases with increasing η and at a high density (η = 0.4) one eventually obtains a bicontinuous phase similar to microemulsions. 54his is apparent from the cluster size distribution, which shows a well defined peak at large s comprising nearly all particles in the simulation box (the inset in Fig. 5d).
Our simulations show that there are no significant changes in the structure within the clusters.This is reflected in the RDF, which exhibits only a slight increase of the first peaks as η increases (Fig. 5c).The hexatic parameter ψ 6 also changes only moderately with η, showing that indeed the structure inside the clusters remains practically unaltered (Fig. 5e).

Percolated network of deposited nanoparticles
Sufficiently far from the critical point, where the decay length of the solvent-induced intersections is smaller than the particle radius, our simulations reveal the formation of fractal-like clusters that consist of thin but elongated aggregates of nanoparticles with nearly crystalline order inside the aggregates (Fig. 3 and Fig. S6 †).At a high nanoparticle density, these aggregates merge and form structures reminiscent of colloidal gels or percolated networks.These structures appear for interaction potentials with deep attractive wells and none or very weak repulsive tails.Increasing the strength of repulsion, which can be done by increasing the Debye length, breaks the percolated network so that the deposited structure consist of separated, rather than inter-connected, but still elongated nanoparticle clusters (Fig. 6).Percolation can thus be induced by manipulating repulsion through varying the salt concentration.
We note that such fractal-like nanoparticle gels are likely only kinetically stable in the sense that, if formed, they can persist for a long time but may eventually collapse into a single compact cluster. 51,55This is unlike cluster fluids or bicontinuous phases, which are thermodynamically stable. 56,57n our simulations, this difference can be seen in the dependence on the initial particle configuration.In the former case, starting from a densely packed, rather than a random configuration, leads to a single crystalline cluster (with defects); in the latter case the densely-packed aggregate breaks and the system evolves towards the structure consistent with the one obtained with the initially random configuration (Fig. S11 †).
The possibility to control the formation of 2D particle networks and tune their percolation properties is exciting and may find applications in optoelectronic devices and sensors. 58

Deposition of quantum dots on a silicon substrate
To test the predictions of our simulations, we performed experiments on deposition of quantum dots (QDs) driven by the interplay of the Debye-screened electrostatic repulsion and the critical Casimir attraction.The setup is shown in Fig. 7a.QDs were synthesized following the reported procedure. 59The native hydrophobic oleate ligands were replaced with hydrophilic 2,3-dimercapto-1-propanesulfonates (DMPS).The diameter of the QD was 2a = 5.5 nm ± 13% as measured from small-angle X-ray scattering.The QDs were dispersed in an offcritical (30 w/w%) water-2,6-lutidine mixture, which was placed in contact with a silicon substrate made hydrophilic by oxygen plasma treatment.The substrate was in contact with a thermostated steel block and heated to temperature T, such that ΔT = T c − T > 0, where T c is the critical temperature of the water-lutidine mixture.Keeping the bulk of the solvent at room temperature drives a vertical temperature gradient, making the surface of the substrate the only location for QD deposition by critical Casimir forces.After one hour the substrate was removed, dried, and imaged with scanning electron microscopy.
We measured the size of solvent density fluctuations, the correlation length ξ, in proximity of the substrate prior to each deposition by measuring the hydrodynamic diameter of the fluctuations with dynamic light scattering and correcting for the T-dependence of the viscosity. 60As the probe beam of the light scattering setup had a diameter of ∼0.5 mm, it measured the bulk correlation length.Thus, the measured ξ describes the decay of the critical Casimir interactions in the bulk fluid, rather than in the immediate nano-scale environment of the substrate, as in the simulations.In the experiments, we consider the values of ξ ranging from 10 to 19 nm (Fig. 7b), which in terms of the QD radius is 3.6a to 7a.Using the results of our MC simulations (Fig. 1h), we estimate that this bulk decay length corresponds to modified decay lengths between 1.5a and 2a at the substrate.
We obtained QD structures that are consistent with the simulation results.Fig. 7c shows micrographs of the obtained structures for three correlation lengths ξ and for the Debye length λ D = 10 nm.We observe the formation of clusters that coarsen and evolve into spinodal-like structures upon increasing ξ, as in the simulations (Fig. 3).The lower micrographs in Fig. 7c show amorphous structures of QDs within the clusters, with some tendency for ordering.Unlike in the simulations, which predict hexagonal ordering, the QDs tend to arrange themselves into a square lattice.This is likely because of their cubic-like shapes. 61e also varied the Debye screening length at fixed ξ.Fig. 7d shows that for a fixed correlation length, the size of the clusters decreases as the Debye screening length increases.This is again in agreement with the simulations (Fig. 4), confirming experimentally the coarsening of the spinodal structures and thus demonstrating the control offered by the critical Casimir interactions over the deposited nanoparticle structures.
While we observe the same trends in the simulations and experiments, including the appearance of liquid-like and spinodal structures of clusters, the simulations additionally predict the formation of fractal-like gels and bicontinuous phases, which we could not identify in the experiments.The absence of nanoparticle gels is not very surprising, as our goal with experiments was to explore the critical regime, while gels appear farther away from the criticality, which gives rise to strong shortranged interactions between the particles (Fig. 6).Bicontinuous phases and long-range cluster fluids appear in a relatively narrow range of parameters and are not easy to locate even with MD simulations.A more detailed investigation is thus desirable, particularly with MC simulations or classical density functional theory, which are more efficient in determining thermo- dynamically stable phases.However, we also cannot exclude that additional effects might be at play in the experiments, which are not taken into account in our model, but which may affect the formation of such mesophases.In particular, the electrostatic interactions were taken in the simplest DLVO form that does not account for solvent's critical fluctuations.It is possible that the varying correlation length influences the nanoparticle surface charge density and hence the strength of electrostatic repulsion. 61QD's polydispersity and ligands may further affect the stability and kinetics of mesophases.These effects require separate experimental and theoretical investigations and fall beyond the scope of the present work.

How interactions define structure
Our multiscale simulations and experiments suggest that there is a relation between the length scales of the interactions and the deposited nanoparticle structures, which we discuss below.The variety of the structures obtained by the simulations are summarized in Fig. 8.
• Increasing λ D leads to a larger energy barrier and a stronger repulsive tail in the interaction potential, which breaks large aggregates and gives rise to smaller dispersed clusters.Examples: λ D /a = 10 in Fig. 4 and 7c.
• Spinodal-like morphologies appear for sufficiently large λ C ≳ a and likely separate the dispersed and coarsened structures.Large λ C and λ D lead to a decreased potential depth (Fig. 2) and hence the spinodal structures are more amorphous than the coarsened ones.Examples: λ C /a = 3 in Fig. 3 and ξ = 19 nm in Fig. 7.This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.
• Cluster fluid or bicontinuous phases (depending on the density) emerge when λ D ≳ λ C ≳ a, provided there is a weak repulsive tail in the interaction potential.An example is shown in Fig. 5.Such mesophases appear in a narrow parameter range, and we could not easily locate them.For increasing λ C , the attraction becomes weaker (Fig. 2) and hence these phases transform into the dispersed ones.For decreasing λ C , the attraction increases, but the long-range behaviour becomes determined by the electrostatic repulsion, leading again to the dispersed structures.
• Fractal-like nanoparticle gels or percolated networks appear for λ C < a and λ D ≲ a, provided the depth of the interparticle potential is sufficiently large so that a nanoparticle attaches to an existing cluster virtually irreversibly.Such strong shortrange potentials exist in the wetting regime and are caused by capillary bridging. 34,37,38In contrast to the critical regime, these forces be well-controlled also for large colloids. 62xamples: λ C /a = 0.53 in Fig. 3 and 6.
In addition to the interplay of the length scales of the interactions, also the potential depth and the height of the repulsive tail (energy barrier) play an important role in the formation of these structures.The strength and the range of critical Casimir interactions cannot be varied independently, but the repulsive interactions can be additionally controlled by changing the nanoparticle charge, which has been assumed the same throughout this work.Clearly, the parameter range, where various structures appear, may be modified by the particle charge.However, we expect the relations between the decay lengths to hold for the structures emerge due to the interplay of the long-range tails of the attractive and repulsive components.
We finally note that the repulsion can also be induced by other means.For instance, for nanoparticles with magnetic cores, the repulsion can be controlled with an external magnetic field, 63,64 which may provide a convenient way to steer the structural transformations without the need to change an experimental setup.

Conclusion
We have developed an experimental method and a multiscale simulation framework to study the behaviour of charged hydrophilic nanoparticles deposited on hydrophilic substrates, immersed in a near-critical binary liquid mixture intermixed with salt.With Monte Carlo simulations of the Ising model, which mimics the binary mixture, we found that the critical Casimir interactions between nanoparticles can be drastically reduced at the substrate (Fig. 1).Notwithstanding, our multiscale molecular dynamics simulations predict the formation of a large variety of superstructures that can be controlled by varying composition of the mixture or temperature and salt concentration (Fig. 3-6 and 8).In particular, the simulations show that the structures coarsen with increasing solvent correlation length and decreasing Debye screening length.The experiments support these trends (Fig. 7), demonstrating new opportunities to steer nanoscale structures with external control.These results signify that nanoparticles in near-critical electrolytic solvents represent an exciting system with the potential for controlled fabrication of complex structures.We hope that our work will stimulate further experimental and theoretical studies that may lead to the emergence of novel fabrication techniques and optoelectronic devices, particularly those based on quantum dots. 65

Multiscale simulations
Interactions between deposited nanoparticles.The interaction potential between deposited nanoparticles has been taken as a sum of three contributions where U rep , U el and U C are the interaction potentials describing the short-range and electrostatic repulsion and critical Casimir attraction, respectively, and r is the center-to-center distance between the particles.Since we focus on nanoparticles with the dielectric constant comparable to that of a critical solvent, we did not include the van der Waals (dispersion) interactions in eqn (2).For non-vanishing dispersion forces, one will have to adapt the electrostatic repulsion to partially compensate for the der Waals attraction by re-adjusting the surface charge or salt concentration.Short-range Lennard-Jones repulsion.For the short-range repulsive interactions we took 66 where A rep = ε LJ π 2 σ 6 /(630a where Z is the particle charge (in units of the elementary charge), l B = βe 2 /(ε r ε) is the Bjerrum length and k D = λ D −1 = (4πρl B ) 1/2 is inverse of the Debye screening length λ D (all in Gaussian units).Here ρ is the salt concentration, ε the vacuum permitivity and ε r the relative permitivity.We took ε r = 7 corresponding to a critical water-lutidine mixture, 42 which gives l B ≈ 8 nm at room temperature.We used the charge Z = 2.13 throughout this work.The values of λ D and A el used in MD simulations are shown in Table S2.† It has been shown that for 2,6-lutidine-water mixture the effect of various salts (viz., potassium chloride, potassium bromide and magnesium nitrate) on the correlation length ξ and the critical temperature T c was minor. 67We correspondingly assumed that adding salt does not change the critical Casimir interactions.
Critical Casimir attraction.We used Monte Carlo simulations of the Ising model to compute the critical Casimir potential U C between the nanoparticles.The Ising model mimics an incompressible binary liquid mixture or a simple fluid.Its Hamiltonian is given by where σ i = ±1 is a classical spin at lattice site i so that σ i = ±1 corresponds to occupation of site i by one of the species of a binary liquid mixture, and J is the coupling constant.The sum 〈ij〉 runs over all neighbouring pairs of spins.The spins belonging to the colloids (and to the substrate, if present) were frozen, so that the summation in eqn ( 6) runs only over nonfrozen lattice sites, corresponding to the space accessible to the solvent.We considered the so-called normal surface universality class, 68,69 which corresponds to infinite surface fields and implies that effectively the surface spins are also fixed during a simulation.Casimir interaction potentials between two spherical particles were computed by using the local field approach. 11,70n all simulations, the system size was 200 × 60 × 60 lattice units in the x, y and z directions, respectively (Fig. 1e and g).The particle diameter was 11 lattice sites.The averaging was performed over 2 × 10 6 hybrid Mont Carlo steps 29 with preliminary thermalization of 2 × 10 5 MC steps.For the configuration with the maximal distance between the particles (89 lattice sites), the number of MC steps was ten times larger to achieve a better accuracy (this value was subtracted from the Casimir potentials obtained at other distances to have zero at large separations).
Molecular dynamics simulations.We wrote our own program to perform 2D molecular dynamics simulations of deposited nanoparticles.In this program, we used the standard Velocity-Verlet integration algorithm.Particles were subject to random thermal forces generated by the solvent.The time step was 10 -4 τ d where τ d = a 2 /D is the diffusion time for the diffusion coefficient D = k B T/γ (a is the particle radius).The friction coefficient γ = 10 was used in all simulations.The averaging was performed over the time 100τ d .The cut-off for the force acting between two particles was 30a to account for the long-range character of the interactions.The initial configuration of nanoparticles was random unless otherwise stated.
Analysis.From the obtained trajectories we determined clusters of nanoparticles.A particle was considered as part of a cluster if the distance to at least one of the cluster's particles was less then 2.5a.Within each cluster, we calculated the order parameters Here N is the total number of nanoparticles, 〈⋯〉 means thermal average, and where n i is the number of neighbours of the ith nanoparticle within its cluster, and θ j,i,j+1 is the angle between nanoparticles i, j and j + 1.We computed ψ 4 and ψ 6 , which describe quadratic and hexagonal order (ψ 4 = 1 and ψ 6 = 1 for perfect square and hexagonal lattices).However, in all cases considered, we have not found any traces of quadratic order (Fig. S12 †).In the main text, therefore, we presented the results only for the hexatic parameter ψ 6 .

Experimental
We synthesized oleate-passivated CdSe quantum dots (QDs) by following procedures reported in the literature. 59After washing, we redispersed the QDs in 10 mL of hexanes, and measured their concentration by following the sizing curve reported in the literature. 71We evaluated the size distribution of the QD dispersion by fitting the small-angle X-ray pattern to a spherical form factor convoluted with a normal distribution centered around the radius a and of width σ.The polydispersity was then extracted as (σ/a) × 100%.The polydispersity can also be evaluated from spectrophotometry by measuring the half-width at half-maximum of the first exciton peak and following the fit provided in the literature. 72e replaced the native hydrophobic oleate ligands with hydrophilic DMPS ligands by following a ligand-exchange procedure.A 40 μM dispersion of QDs in hexane was placed in contact with a 40 mM solution of DMPS in N-methylformamide (NMF) and stirred vigorously for 1 hour.During this time, the QDs migrated from the apolar to the polar phase, as visually confirmed by the apolar phase becoming colorless and the polar phase gaining color.We then stopped the stirring, allowed for the two phases to separate, and transferred the polar phase to a centrifuge tube.Then we precipitated the QDs by adding a volume of acetonitrile equal to the volume of NMF and mixing.After centrifuging at 3000g for 20 minutes, we discarded the clear supernatant and redispersed the QDs at a concentration of 2 μM in the binary solvent, a 30% w/w solution of 2,6-lutidine in water, with an additional 0.89 mM sodium chloride to control the Debye screening length and 10 mM sodium hydroxide to improve colloidal stability of QDs.
We transferred the colloidal dispersion to a glass vessel for dynamic light scattering, where we placed it in contact with an oxygen plasma-treated silicon substrate placed at the top and attached to a steel block.Placing the substrate at the top allows us to rule out sedimentation as a contribution to the structures we observe and minimizes convection.The substrate was heated by means of a heating element embedded in the steel block, to a temperature T, such that ΔT = T c − T > 0, where T c is the critical temperature of the binary solvent.By keeping the bulk of the solvent at room temperature (RT), we drive a steep vertical temperature gradient going from T in proximity of the substrate to RT far from the substrate (Fig. 7a).To ensure experiment reproducibility, we measured the size of the solvent density fluctuations prior to each deposition in proximity of the substrate by measuring the hydrodynamic diameter of the fluctuations with dynamic light scattering and correcting for the T-dependence of the viscosity. 60he temperature range was 20 °C to 33.1 °C.After 1 hour, the substrate was retrieved, dipped in pure acetonitrile, and dried under vacuum prior to electron microscopy studies.By controlling the Debye screening length through the concentration of salt dissolved in the binary solvent, and the solvent correlation length through temperature, we are able to drive the assembly of nanoparticles towards structures of strikingly different morphologies. 61

Fig. 1
Fig. 1 Critical Casimir interactions between nanoparticles.(a) Schematics of two nanoparticles at a substrate immersed in a binary liquid mixture.The nanoparticles and the substrate favour the same component of the mixture.(b and c) Schematic bulk phase diagrams of a binary liquid mixture, around its lower critical point, and of the Ising model (eqn (6) in Methods), which models the binary mixture close to criticality.Concentration of the binary mixture corresponds to magnetization in the Ising model and the difference in the chemical potentials of the two components of the mixture to bulk field h.The thin vertical and horizontal lines show the paths used in our experiments (b) and simulations (c).(d) Critical Casimir potential for two nanoparticles in bulk and at the substrate obtained by Monte Carlo (MC) simulations of the Ising model for the values of bulk field h and reduced temperature τ = (T − T c )/T c indicated on the plot.The shaded area shows the region of steric exclusion.Both particles are located at z/a ≈ 1.45 above the substrate (unless they are in bulk).(e) Magnetization heat maps from the MC simulations.Nanoparticle separation r/a ≈ 4.91.(f and g) The same as (d and e) but for h and τ further away from the criticality along the path of constant magnetization (m = 0.2, cf.panel (i)).The white lines mark zero magnetization.(h) Decay length λ C of the critical Casimir potential for two nanoparticles at the substrate as a function of λ C in bulk along the paths of constant magnetization m = 0.2 and m = 0 (critical composition).(i) Thermodynamic paths defined by constant reduced temperature τ, constant bulk field h, and constant bulk magnetisation m.The colored symbols show the values of the decay length λ C of the critical Casimir potential for two nanoparticles at the substrate.Symbols + and × in (h) and (i) show the values of τ and h used in (d-g).

Fig. 2
Fig. 2 Interactions between nanoparticles at substrate.(a) Example of the critical Casimir, electrostatic and total interaction potentials between two identical nanoparticles at a substrate immersed in a near-critical fluid.The substrate and the particles favour the same fluid component.(b) Diagram showing the parameter space of the Debye screening length λ D and the decay length of the Casimir interactions, λ C .The colored symbols show the values of the depth of the total interaction potential.(c) The same as (b) but showing the energy barrier of the total interaction potential.The diagrams are built along the thermodynamic paths indicated in Fig. 1i.The × symbols in (b and c) denote the values of λ D and λ C corresponding to panel (a).The nanoparticle structures for the paths of constant λ D /a = 1 and constant λ C /a = 2.5 (bulk magnetization m = 0.2) are shown in Fig. 3 and 4 as denoted on panels (b) and (c).The paths of constant bulk field h = 0.007 and reduced temperature τ = 0.007 (at fixed Debye lengths) are shown in Fig. S4 and S5.†

Fig. 3
Fig. 3 Effect of the decay length of the critical Casimir interactions.(a) Total interaction potential between nanoparticles at a substrate immersed in a near-critical fluid for three values of the decay length λ C along the path of constant magnetization m = 0.2 (Fig. 1).The substrate has the same preference for the fluid phase as the nanoparticles.(b) Radial distribution functions (RDFs) and (c) cluster size distributions in systems with the interaction potentials from panel (a) obtained by 2D molecular dynamics (MD) simulations.(d) Average cluster size and hexatic parameter ψ 6 , eqn (7), as functions of λ C .(e) Snapshots from the MD simulations, showing fractal-like, coarsened and spinodal-like structures (see also Video †).In all plots, the Debye screening length λ D /a = 1 and surface packing fraction η = 0.3.See Fig. S3 † for the case of bulk magnetization m = 0 corresponding to the critical composition.For the dependence on η see Fig. S6 and S7.†

Fig. 4
Fig. 4 Effect of the Debye screening length.(a) Total interaction potential between nanoparticles at a substrate in a near-critical fluid for a few values of the Debye screening length λ D measured in units of the particle radius a.The substrate has the same preference for the fluid phase as the particles.(b) Radial distribution functions and (c) distribution of cluster sizes in the systems with the interaction potentials from panel (a) obtained by 2D molecular dynamics (MD) simulations.(d) Average cluster size and hexatic parameter ψ 6 , eqn (7), as functions of λ D .(e) Snapshots from MD simulations, showing coarsened, spinodal-like and dispersed structures.In all plots the reduced temperature τ = 0.007 and the bulk field of the Ising model h = 0.003 (giving λ C /a ≈ 2.5), and the surface packing fraction η = 0.3.For the dependence on η see Fig. S8-S10.†

Fig. 5
Fig. 5 Effect of nanoparticle density.(a) Electrostatic, critical Casimir and the total interaction potential between two nanoparticles at a substrate in a near-critical fluid.The substrate and the particles prefer the same fluid component.(b) Snapshots from MD simulations for a few values of the two-dimensional packing fraction η = πNa 2 /A, where N is the number of particles, A the surface area and a the particle radius.(c) Radial distribution functions and (d) cluster size distributions for the same systems.The cluster size distributions are averages of 100 simulations.(e) Average cluster size and hexatic parameter ψ 6 , eqn (7), as functions of η.The reduced temperature τ = 0.007 and bulk field h = 0.007, giving λ C /a ≈ 2, and the Debye screening length λ D /a = 4.

Fig. 6
Fig. 6 Percolated network of deposited nanoparticles.(a) Interaction potential between deposited nanoparticles in a near-critical fluid for two values of the Debye screening length λ D .The substrate and the particles prefer the same fluid component.(b) Average cluster size and hexatic parameter ψ 6 , eqn (7), as functions of λ D .The vertical dash line shows approximately the 'percolation' value of λ D .(c) Snapshots from MD simulations for the potentials from panel (a).The reduced temperature τ = 0.1 and the bulk field h = 0.01 of the Ising model give the bulk magnetization m = 0.2 and the decay length of the critical Casimir interactions between the deposited nanoparticles λ C /a ≈ 0.53.

Fig. 7
Fig. 7 Deposition of quantum dots on a silicon substrate.(a) Schematics of the experimental setup.Hydrophilic quantum dots (QDs) of diameter 2a = 5.5 nm are immersed in a water-2,6-lutidine mixture at off-critical composition (30 w/w%) in contact with a hydrophilic silicon substrate.Temperature varies from room temperature in the bulk to a nearly critical at the substrate.(b) Diagram in the plane of the Debye screening length λ D and bulk correlation length ξ, showing the values of λ D and ξ used in the experiments (circles).The bulk correlation lengths correspond to the decay lengths of the Casimir interactions between QDs at the substrate ranging from about 1.5a to 2a (Fig.1h).The squares indicate the simulations results shown in Fig.3 and 4. (c) SEM images of deposited QDs for a few values of the bulk correlation length ξ and for the Debye screening length λ D = 10 nm.The deposited structure coarsens with increasing ξ, similarly as in the simulations (Fig.3).(d) SEM images for two values of the Debye screening length and for the bulk correlation length ξ = 10 nm.The size of clusters decreases with increasing λ D , in agreement with the simulations (Fig.4).In (c and d), the lower row shows the enlargement of small areas from the upper row, demonstrating crystal-like structures inside the clusters.The sizes of bars are 125 nm (upper row) and 25 nm (lower row).

Fig. 8
Fig. 8 Summary of nanoparticle structures.Diagram shows the summary of structures in the plane of the Debye screening length λ D and the decay length of the critical Casimir interactions between nanoparticles of radius a deposited on a substrate.Examples of the coarsened and dispersed structures are shown in Fig. 4e, λ D /a = 0.5 and λ D /a = 10, respectively.Spinodal structures are shown in Fig. 3e (λ C /a = 3) and Fig. 4e (λ D /a = 2.5).Cluster fluid and bicontinuous phases are presented in Fig. 5. Examples of percolated and non-percolated nanoparticle networks are shown in Fig. 6.
6) and x = r/a, with a being the particle radius, and ε LJ and σ the parameters of the intermolecular Lennard-Jones potential.Taking a = 2.75 nm as in the experiments, and using typical values ε LJ = 0.4 kcal mol −1 and σ = 0.3 nm, gives A rep ≈ 1.1 × 10 -8 kcal mol −1 ≈ 1.78 × 10 -8 k B T at room temperature, which has been used in all simulations.Electrostatic repulsion.The electrostatic interactions were taken in the DLVO form 1 βU el ðxÞ ¼ A el e ÀakDðxÀ2Þx ;