Temperature-induced structural transitions in self-assembling magnetic nanocolloids †

With the help of a unique combination of density functional theory and computer simulations, we discover two possible scenarios, depending on concentration, for the hierarchical self-assembly of magnetic nanoparticles on cooling. We show that typically considered low temperature clusters, i.e. defect-free chains and rings, merge into more complex branched structures through only three types of defects: four-way X junctions, three-way Y junctions and two-way Z junctions. Our accurate calculations reveal the predominance of weakly magnetically responsive rings cross-linked by X defects at the lowest temperatures. We thus provide a strategy to ﬁne-tune magnetic and thermodynamic responses of magnetic nanocolloids to be used in medical and microﬂuidics applications.


Introduction
The first stable suspensions of magnetic nanoparticles, nowadays known as ferrofluids 1,2 , were synthesised more than 50 years ago 3 . Since then, the design of magnetic colloids has steadily progressed: it is now possible to functionalise the nanoparticle surfaces 4,5 , tune their shape [6][7][8][9][10][11][12] and control their internal crystalline structure [13][14][15] . Such a flexibility in design has opened broad perspectives for applications of magnetic nanocolloids in biomedicine, such as sensors, contrast agents, drug targeting, hypothermia cancer-treatment therapy, cell manipulation and many others [16][17][18][19][20][21] ; they also improved various technologies, such as loudspeakers, seals, dampers, purification and separation devices [22][23][24] . Magnetic nanoparticles are attractive for all aforementioned applications because of their single-domain nature. Indeed, magnetic nanoparticles have intrinsic magnetic moments and behave as nanomagnets as long as their size does not exceed the size of a magnetic domain. The latter depends on the ferri-or ferromagnetic material the particle is made of, but it is usually in the range of 10-50 nm. Another important common property of all these magnetic nanocolloids is that, once in solution, they undergo brownian motion. Moreover, if the magnetic crystallographic a University of Vienna, Sensengasse  anisotropy of a nanoparticle is energetically strong enough, then the orientation of the magnetic moment is coupled to the orientation of the particle as a rigid body 25 , so that the two cannot rotate independently.
Notwithstanding the way the surface of a magnetic nanoparticle is treated (be it special funtionalisation, steric or electrostatic stabilisation) there is an inevitable long-range, anisotropic dipolar interaction between nanoparticles' magnetic moments. This interaction introduces a preferred headto-tail orientation of dipole moments and, if strong enough to compete with thermal fluctuations, leads to a directional self-assembly of particles in linear flexible chains [26][27][28][29][30][31] . Of course, controlling self-assembly in magnetic nanocolloids is not restricted by magnetic forces only. Van der Waals, electrostatic or chemical interactions of functionalised particle surfaces, solvophilicity/solvophobicity of nanoparticle shells; all these and many others can result in versatile selfassembly scenarios [32][33][34] . It is interesting to observe in this respect the analogies between the anisotropically interacting dipolar nanoparticles and patchy colloids with a reduced valence 35,36 . Indeed, restricting the particle-particle relative orientations that result in a bond, leads to self-assembly processes that resemble the one occurring in dipolar nanoparticles, suggesting the possibility to extend the recent investigations of patchy colloids to dipolar fluids. For example, many novel and technologically-relevant phases, such as open lattices 37,38 , equilibrium gels 39,40 and reentrant gas-liquid phase separations 41,42 , can be produced by a fine tuning of the patchy colloids parameters.
Along with an impressively rich spectrum of phases and structures induced by the aforementioned short-ranged made of two rings cross-linked by one X defect. Fig. 2 provides a cartoon of these different aggregation pathways.

Model
We investigate systems composed of N magnetic nanoparticles at varying temperature T and concentration c. Two particles i and j interact through the potential V tot = V HS + V dd , where V HS is the hard-sphere potential given by where r i j = | r i j | is the inter-particle distance and d is the particle diameter. V dd is the dipole-dipole interaction, defined as follows: where µ k is the magnetic moment of particle k;ˆ r i j = r i j /| r i j |. All the particle magnetic moments have the same magnitude µ.
In the text, k B is the Boltzmann constant, β = 1/T , lengths are measured in units of the particles diameter d and energy in units of µ 2 /d 3 . Therefore, temperature is measured in units of k B d 3 /µ 2 .

Theoretical Approach
Our detailed analysis of branching in Ref. 56 revealed the possibility to split all possible branching points in the system into three main classes only (see Fig. 1), using the number of ways out w (in the following "valency") as a label. Looking at Fig. 1 one can see that particles of type Z have w = 2, particles of type Y are characterised by w = 3 and particles of type X have w = 4. Each of these defect particles also has a "charge" due to the possible dipolar orientations of defect neighbours: Y ones have charge 1 or -1, since if two dipoles enter the Y particle, then one necessarily exits (charge -1); conversely, if two dipoles exit, one enters (charge +1). By contrast, the charge of Z and X particles is always zero, since if one or two dipoles enter, one or two others always exit from the defect particle.
In our model, defect particles do not interact with each other; rather they interact with "regular" particles. The latter can form chains and rings and as such connect defect particles. One can also think of our system as being composed by four different types of particles: one type can only form chains and rings, the other type (Z) can attach to a ring or to a chain, but cannot connect any of the two to some other cluster, whereas the other two (Y and X) connect chains and rings, effectively serving as cross-linkers with a fixed valency.
In Fig. 1 we showed the nine classes of structures that contain up to one defect particle. The free energy (F) per unit volume of this system can be obtained by summing up the free energies of these classes. In the following, we denote their equilibrium concentrations as g(·).
Here, k B T is the thermal energy. Unlike the monodisperse case 51 , the summation here cannot be carried over the number of particles in a cluster, but one has to introduce a vector i containing the information about not only the size of the cluster but also about its topology, i.e.
where n ≥ 0 is the amount of regular particles in chain segments; (k j ≥ 4) or (k j = 0), j = 1, 2 are the numbers of regular particles in ring segments; m Y , m X and m Z are the numbers of type Y , X and Z defects present in the cluster (in the current study ∑ l=X,Y,Z m l = 0 for primary structures and 1 for the branches). In addition, the partition function of each class also depends on i: Q( i). Finally, v( i) is a characteristic volume, which allows us to express the partition functions in a simple factorised way (see, Table 1). The parameter K( i) is a combinatorial factor, representing the number of entropically distinguishable clusters from the same class, having the same i and with the same partition function Q( i). The free energy can thus be computed by minimising the functional (3) taking into account the mass-balance condition: This condition constraints the total number of particles in a unit volume (φ /v, φ being regular particle volume fraction, v standing for particle volume). Note that the total number of regular particles in a cluster of class i is given by where a defect particle depending on its type contains s Y , s X or s Z regular particles. These parameters are used to simplify the description of the cluster formation. As shown in Ref. 56 , the number of particles in the defects can vary from one to four, however, the actual internal structure is irrelevant. Thus, we can unify the treatment and limit the number of defect particle types. Minimisation using Lagrange multipliers method leads to a solution with a "traditional" functional form: where p is the Lagrange multiplier and µ has the meaning of a chemical potential. The main complexity of the problem is condensed in the analytical calculation of Q( i), due to the long-range magnetic dipole-dipole interactions between nanoparticle magnetic moments. Here, we put forward an approach that makes use of the following assumptions: • each bond is described by an effective free energy: ε L for a linear chain segment, ε R for a ring segment; ε Y L(Y R) for a bond between the member of a chain (ring) segment and a defect particle Y ; ε XL(XR) for a bond between the member of a chain (ring) segment and a defect particle X; ε ZL(ZR) for a bond between the member of a chain (ring) segment and a defect particle Z; ε l , l = X,Y, Z is the free energy loss of forming a defect l compared to a bond between regular particles.
• the partition function Q( i) can be factorised as a product of the aforementioned effective bond free energies and thus be computed as a function of the number of monomers and temperature, for details see Ref 51 ; • the factorisation of Q( i) has a peculiarity, namely, even though it is presented as a product of effective energies, these energies actually do depend on the lengths of corresponding segments; in this way we can consider the interactions beyond the limit of the nearest neighbours; • entropic contributions are calculated with respect to the entropy of chain, which is why the factor 1/k 3ν+1 is included in partition functions of ring containing clusters to capture the difference in entropy between chains and rings arising from the k ways of opening a ring to form a chain, where ν = 0.588 is the self-avoiding random walk exponent; the difference number of self-avoiding paths of chains and rings is proportional to k 3ν 57 .
Building upon past results 51 , we know that the bonding free energies associated to rings and chains are respectively given by where with ζ (3) denoting the Riemann zeta function of three; R n−1 2 stands for the residual of division, and [·] has the meaning of the integer part of the expression in the brackets. The low-T dimer partition function q (note that C(1) = 0, C(2) = 1), derived first by de Gennes and Pincus 58 , is As a first step, we assume that In addition, we introduce three parameters: In the general case x, y, z are slowly changing functions of T * , and from previous investigations 56 we know that x < y < z ≪ 1 for any temperature. Here, however, we neglect their temperature dependence and fit their values by matching the internal energy in simulations and theory, obtaining x = 0.003, y = 0.005 and z = 0.006. We performed an extensive analysis changing the values of these parameters by up to 60 per cent and discovered, that the microscopic characteristics as well as the thermodynamic quantities are not sensitive to such variations of x, y and z.
For the sake of simplicity, we fix s X = s Z = s Y = 1. To justify the latter, we have carefully checked simulation results, finding that the amount of particles in defects is very small in comparison to the number of particles in ring and chain segments. We provide the expression for all the terms in the free energy functional (3) under these assumptions in Table 1.
Using this table, one can minimise the free energy of the system and obtain the Lagrange multiplier p as a function of temperature and nanoparticle concentration.
It is worth noting, that using the approach described above, we can introduce more complicated structures, that can be hierarchically obtained by combining the basic nine elements, as shown in Fig. 3. One can keep adding up defects following this approach to eventually build a percolating network of branched structures. The table can be relatively easily extended for this case, however, the calculation of the K values becomes more cumbersome. In any case, before dealing with the network itself it is crucially important to understand how the initial stage of branching occurs, which structures are more probable at low temperatures and what are the main tendencies in branching for different particle concentrations. In other words, we aim at understanding the second step (the first being the formation of primary structures, namely chains and rings) of the hierarchical self-assembly, thus elucidating the thermodynamics of the aggregation of primarily formed chains and rings into small branched clusters. Physical Chemistry Chemical Physics

Physical Chemistry Chemical Physics Accepted Manuscript
All possible 32 classes (from 10 to 41) of aggregates containing two defect particle. Defects Y are magenta, defects X are blue, defects Z are green, linear segments are black, ring segments are cyan.

Computer Simulations
We carry out Monte Carlo simulations in the canonical ensemble of N = 5000 dipolar hard spheres at different temperatures and concentrations ϕ = N/V ≪ 1, where V is the volume of the simulation box. We simulate down to T * = 0.125 by implementing ad-hoc biased Monte Carlo moves 44,59 . For further details on the simulation approach, see Ref. 44 . We partition particles into clusters by employing a mixed distance/energy criterion: two particles are bonded if their pairinteraction energy is negative and their relative distance is smaller than r b = 1.3 44 , which corresponds to the first minimum of the radial distribution function. Clusters are then classified according to these definitions: • a structure containing two single-bonded particles con-nected by particles having two neighbours is labelled as a chain; • a ring is a cluster containing only particles having two neighbours; • any other kind of aggregate is a branched cluster.

Results and Discussion
The main advantage of the combined approach proposed here is that we can not only scrutinise the topologies of selfassembled magnetic nanoparticle clusters and their distributions at various concentrations and temperatures, but we can also directly follow the hierarchical self-assembly in these systems and estimate the thermodynamic contribution of each cluster type. In order to obtain a clear understanding of the concentration and temperature influence on the self-assembly in magnetic nanocolloids it is convenient to look at the fractions of particles in various clusters. In the following, we use dimensionless units for temperature (T * ) and nanoparticle concentration (ρ * ), as defined in Model and Methods.  decreases, the chains increase in length and eventually the energetic gain of closing into rings becomes large enough for a structural transition to occur (Fig. 4 upper left). Further cooling the system drives the formation of more complex structures, even though their fraction never exceeds that of rings. A different scenario takes place at large nanoparticle concentration. In this case ( Fig. 4 upper right, lower left and right), the environment becomes too crowded for ideal rings to dominate on cooling, and defect structures become increasingly relevant. The temperature of this structural transition grows slowly with growing nanoparticle concentration. We stress that our assumption that clusters can contain up to one defect is valid in a rather broad range of concentrations, as demonstrated by the good agreement between theoretical predictions and simulation results. However, for the highest concentration considered here (Fig. 4 lower right), the one-defect limitation becomes too strong since multiple-defect structures start to appear and the agreement between theoretical and numerical results degrades. It means that the structures considered here are the precursors to the next hierarchical level of self-assembly in magnetic nanocolloids, namely the aggregation of single-defect clusters into complex networks.

Cluster Distribution
Once it is clear that the high-T defect-free clusters are replaced by more complex structures upon cooling, the composition of these defect structures became of primary importance. Our theoretical approach allows us to scrutinise all possible transitions in detail. In order to illustrate the competition between different defect types, we show the temperature and concentration dependence of the fractions of particles in the various clusters in Fig. 5. The colour-code is the same for each plot and is provided on the right: the brighter the colour is, the higher the fraction. The classes shown along the main diagonal (I, II, VII) do not exhibit any maxima in the range of studied parameters, whereas all the other structures have high fractions in well-defined T * − ρ * regions. In this figure one can clearly follow the hierarchical nature of the self-assembly in the systems under study. For example, the intracluster defects (VIII, IX) start to emerge with decreasing T * , but then quickly get replaced by more complex aggregates (III, IV, V, VI). Note that all these transitions are caused by the magnetic dipole-dipole interaction and entropy competition, which becomes more and more intricate.
Y structures first start to emerge and then slowly disappear in favour of X-structures. At the same time, rings with internal defects and tennis-racket-like structures (IV) are clearly overtaken by rings with two chain tails. If we cool the system further down, basically all defect-cluster fractions exhibit a clear maximum and start rapidly decreasing. For low concentrations, the majority of magnetic nanoparticles is aggregated in rings. For higher concentrations, though, the fraction of particles in rings has a maximum, signalling a transition towards the formation of higher-order defect clusters such as those in the VII class, i.e. double rings. These double rings thus form at temperatures lower than single rings, but at relatively high nanoparticle concentration they clearly become the dominant class. Another interesting aspect shown by Fig. 5 is that the fractions of chain-only defects decrease faster on cooling than those with at least one ring segment. This is a clear reminiscence of the low-concentration structural transitions on cooling, where the complete redistribution from chains into rings is observed in the same temperature range 51 .

Thermodynamics
Next we investigate the effect of the observed structural transitions on the thermodynamics of the magnetic nanoparticle systems. We start by considering how temperature and particle concentration affect the internal energy per particle, shown in Fig. 6. Several conclusions can be drawn when looking at Fig. 6. At any concentration the internal energy has an inflection point, which generates a maximum in the T dependence of the constant volume specific heat C V 60 . The position of the maximum in C V , as shown in the inset, shifts towards low temperatures as the concentration decreases. Independently from concentration, there is only one specific heat maximum, whose position is in close correlation with the temperature of the first structural transition, be that chain-ring or chain-defect transition (compare to Fig. 4).
Taking advantage of the theoretical approach, one can compute the free energy and estimate its change across the various structural transitions. First we estimate the thermodynamic driving force which leads to the ring formation by calculating the free energy difference between a system composed by only chains and a system with chains and rings. At low concentrations, the possibility of forming rings drives to a dramatic decrease of the free energy and the contribution of rings to the system free energy becomes dominant at low temperature. These results are plotted with dashed lines in Fig. 7. For denser systems, though, the free energy gain due to ring formation is marginal. As a second step, we consider the free energy difference between a system composed only by chains and rings and a system composed of chains, rings and defects. The resulting driving force for branching is plotted with solid Theoretical solid lines agree well with simulation results, whereas the dashed curves, which show the internal energy of the system in which only chains can form, strongly deviate. The internal energy has a prefactor, which is set to unity in this plot: vd 3 /µ 2 V , where v and V are the volumes of the particle and the system respectively, µ is the value of the nanoparticle (diameter d) dipole moment. In the inset the temperature of the specific heat maxima T m , i.e. the positions of the internal energy inflection points, is plotted as a function of the nanoparticle concentration. lines in the same figure (Fig. 7). The results show that, indeed, there is a significant driving force for branching, but only at intermediate and high densities in the studied range of temperature. This result confirms that the loss in the translational entropy of primary clusters on branching becomes less-and-less relevant with increasing density in comparison to the energetic contribution of the additional interaction.

Conclusions
In order to control the self-assembly of magnetic nanoparticles, one needs to gain a fundamental understanding of the interparticle interaction. Here, we reveal two scenarios of hierarchical self-assembly induced by cooling. Magnetic dipoledipole interactions leads to the formation of non-bulky lowdensity structures in which the head-to-tail or antiparallel orientations of the dipoles are favourable. This can be realised by forming chains, rings, or various branches. The interplay between energy gains and entropy losses determines the internal structure of the self-assembling system and, as such, drastically influences the macroscopic responses of the systems. Even though nanoparticle rings have lower energy, at moderate temperatures, linear chains are predominant owing to their higher entropy. This changes on cooling. As a result, Fig. 7 Free energy difference in units of T * calculated per particle. Dashed lines describe the free energy gain due to the ring formation in the model system with chains and rings only 51 . The free energy gain obtained through the formation of branched structures (within the framework of the model proposed here) is plotted with solid lines. The legend for different nanoparticle concentrations is provided.
when the system of magnetic nanoparticles is highly diluted, the rings with a vanishing dipole moment turn out to be dominant. The same happens to the Y-type junctions, which form in the systems at higher nanoparticle concentration. However, the entropy gain of forming a Y-junction does not compensate for the energy decrease gained from the ring formation, thus making the structural behaviour of highly magnetic nanocolloids much more complex than expected. In this work, using a combined novel analytical-numerical approach, we discover a series of fascinating structural transitions. Indeed, we find that the system of magnetic nanocolloids exhibits two different scenarios of the hierarchical self-assembly. At low particle concentrations, first highly magnetically responsive linear chains are formed, which further close in rings, thus making the whole system barely susceptible to an external magnetic field. The weak interaction between nanoparticle rings increases on cooling, which leads to the formation of a particular type of defects, playing the role of ring cross linkers. In contrast, at higher nanoparticle concentration, at moderate temperature, chains are replaced by entropy-induced branched structures, which can be described by only three types of defects, X, Y and Z, classified according to the number of linear segments connected by the defect. On cooling Y and Z defects disappear, and the X-type defects dominate, merging together nanoparticle rings, keeping the system weakly magnetically responsive. The defect structures observed here can be considered as the precursors to the network formation, ushering in the next step of the hierarchical self-assembly of magnetic nanocolloids. Our theoretical frame-work developed here, combined with the approach put forward in Ref. 54 can serve as a starting point to the description of the network formation.