Thermophoretic torque in colloidal particles with mass asymmetry

We investigate the response of anisotropic colloids suspended in a fluid under a thermal field. Using nonequilibrium molecular dynamics computer simulations and nonequilibrium thermodynamics theory, we show that an anisotropic mass distribution inside the colloid rectifies the rotational Brownian motion and the colloids experience transient torques that orient the colloid along the direction of the thermal field. This physical effect gives rise to distinctive changes in the dependence of the Soret coefficient with colloid mass, which features a maximum, unlike the monotonic increase of the thermophoretic force with mass observed in homogeneous colloids.


I. INTRODUCTION
Thermal gradients induce mass and charge coupling effects that can potentially be employed to construct devices for energy conversion applications. Since the early observations of thermodiffusion by Ludwig and Soret [1,2] in the 19th century, other coupling effects have been reported over the years. Lehmann showed, shortly after the discovery of liquid crystals (LCs), that cholesteric LCs adopt a uniform rotation as a response to thermal gradients [3]. Peltier and Seebeck demonstrated the coupling between electric currents and thermal gradients [4], which is the basis of thermoelectrics [5]. More recently, it was shown that thermal gradients induce the polarization of liquid water (thermal polarization, TP) and thermal orientation (TO) in molecular fluids of anisotropic particles [6][7][8]. These works highlighted the importance of particle anisotropy in molecules and colloids as an important variable driving thermal orientation effects.
The response of axially symmetric particles, e.g., spherocylinders, to thermal gradients was investigated in Refs. [9,10]. These works used kinetic theory to study particles suspended in a gas (aerosols) under thermal gradients. It was found that the thermophoretic drift is anisotropic. This observation is relevant in anisotropic particles [11], since the anistropy could influence deposition processes in suspensions. More recently, the thermophoretic drift of spherocylinders has been investigated using hydrodynamic computer simulations [12]. These authors reported anisotropic thermophoresis in fluid suspensions as well, and concluded that the anisotropy does not induce particle orientation in a homogeneous thermal field. We discuss in this work how inhomogeneous mass distribution inside a colloid can impart an orientation in such homogeneous fields.
The manipulation of suspensions of anisotropic particles with thermal gradients is attracting significant attention and motivating new experiments. The thermodiffusive behavior of * j.olarte@imperial.ac.uk † mrubi@ub.edu ‡ f.bresme@imperial.ac.uk colloidal rods (fd viruses) was studied recently [13], and it has been shown that the heating of metallic nanorods in a polarized optical trap induces a thermal gradient around the rod and a torque that can be significant (∼10 2 pN nm) [14]. Theoretical and experimental studies of anisotropic particles under thermal gradients have uncovered fascinating behaviors that could open new routes to manipulate small colloids using thermal fields. Understanding the behavior of these particles under nonequilibrium conditions is of fundamental interest to explain and predict coupling effects out of equilibrium.
In this work we investigate the response of small anisotropic colloids dispersed in a dense fluid, which is subjected to a thermal gradient. We will analyze the thermal coupling effects arising from the colloid anisotropy, when the colloids feature an inhomogeneous mass distribution. We will demonstrate that mass anisotropy, which is a general feature of molecular assemblies and a property that can be tuned in colloids with heterogeneous composition, does couple with the imposed thermal fields. The coupling leads to a transient thermophoretic torque and the colloids adopt in the stationary state a preferred orientation. We also show that this coupling influences and modulates the Soret coefficient of the anisotropic colloids.
The systems we are interested in, colloidal suspensions of anisotropic particles with asymmetry in the mass distribution, involve high fluid densities and short mean free paths, the latter being smaller than the characteristic particle size. Further, frictional effects, hydrodynamic interactions, thermal fluctuations, and the fluid structuring around the particle become important in describing the thermophoretic behavior of the suspension. These issues limit the applicability of kinetic theory. We therefore used nonequilibrium atomistic molecular dynamics to perform our investigation. We also use nonequilibrium thermodynamics theory [4,8] to derive phenomenological equations that describe the physical effects discussed herein.

II. METHODOLOGY
We performed simulations of a single colloid suspended in an atomic fluid at a density characteristic of a liquid. The fluid was modeled using the WCA potential [15] at an average density of ρ = (N/V )σ 3 = 1.0, where N is the number of solvent particles, V = L x × L y × L z is the volume of the simulation cell, and σ is the diameter of the solvent particles. The colloid was modeled as a rigid chain of tangent spheres using a collection of beads of diameter, σ . We considered both symmetric and asymmetric mass distributions. The latter was achieved by changing the mass of one bead at the end of the colloid, such that the total mass is given by m total = (N − 1)m + m end . The mass, m, of some of the beads is equal to the mass of the solvent particles, m s , and it is defined as m = m s = 1 hereafter. For the symmetric case we used m total = (N − 2)m + 2 × m end . The degree of mass asymmetry can be quantified by introducing a new quantity, the where r g is the geometric center of the colloid, m i and r i are the mass and position of the beads in the colloid along the colloid axis, and N c is the number of beads in the colloid. For the model considered in this work the dipole is defined as with γ = m end /m and γ = 1 (d = 0) for the asymmetric and symmetric cases, respectively. We show in Fig. 1 a typical simulation cell illustrating the position of the thermostats employed to set up the thermal gradient (see also Supplemental Material [16]). The center of mass of the colloids was tethered using a harmonic potential (see below) to the geometric center, r c , of the left and right containers in the simulation cell (see Fig. 1), while the colloids rotate freely around their center of mass. Colloids of different lengths and different values of the γ parameter defined above are shown in Fig. 1.
The WCA potential was used to model the interparticle interactions, where represents the interaction strength and σ is the diameter of the solvent particles. The interaction parameters were = σ = 1 for all colloid-fluid and fluid-fluid interactions, while the atoms within the nanoparticle did not interact with each other. We use σ and to define the usual Lennard-Jones units. 16000 particles were initially placed in a simulation box of size (40,20,20)σ , for an average reduced density ρ = 1.0 and reduced temperature T = 2.5. Two colloidal particles were placed at positions (10,10,10)σ and (30,10,10)σ and random orientations. Each colloid was treated as a rigid body. The harmonic potential to restrain the translational motion is defined as V (r) = k(r − r c ) 2 , where k is the force constant and r c is the position of the reference point at the beginning of the simulation, corresponding to the geometric center between the hot and cold thermostats (see Fig. 1). The force constant was set to kσ 2 / = 10 3 . A thermal gradient and a heat flux were imposed on the system by using boundary driven nonequilibrium molecular dynamics (BD-NEMD) simulations, using a time step of δt = 0.0025. All the trajectories were integrated using LAMMPS [17]. Thermostatting regions with dimension (4,10,10)σ were set up at the edges and in the center of the prismatic box (see Fig. 1). The motion of the fluid atoms was integrated using the velocity Verlet method, and the translational and rotational degrees of freedom of the colloids were integrated using the method of quaternions for rigid bodies. The system was equilibrated at the desired average temperature of the NEMD run by using the Nosé-Hoover thermostat. We then applied the thermal gradient by thermostatting the hot and cold regions at the desired temperatures. The evolution of the temperature gradient through the system was monitored periodically over a simulation of 10 6 steps. The gradient reaches the stationary profile before this period. Production simulations were performed under the same conditions for an additional 10 8 time steps and the orientation vector, n, of each colloid was sampled every 100th step.
The thermophoretic force on each colloid is computed as a function of the thermal gradient, by sampling the average displacement of the center of mass of the colloid and multiplying by the force constant of the applied harmonic potential. The total torque is referred to the center of mass of the colloid, computed for symmetric and asymmetric colloids as a function of the orientation of the colloid, but without integrating the rotational degrees of freedom. The systems were equilibrated for 10 6 steps, and the torque was averaged every 100 steps for 10 7 steps and 30 different runs.

III. RESULTS
Thermal gradients induce a thermophoretic force (TF), F T , on the colloids, and the latter drift either to the hot (thermophilic) or cold (thermophobic) regions. The TF can be calculated by monitoring the force on the colloids exerted by the solvent. The TF can be related to the Soret coefficient, S T , through [18,19], F s − F T = S T k B T ∇T , where T is the local temperature, ∇T is the thermal gradient, and F s and F T are the thermophoretic forces on the solvent and colloid, respectively. This approach has been successfully tested before to compute the Soret coefficient of colloids and binary mixtures [18,19]. In our case the thermophoretic force dominates the value of the Soret coefficient. The contribution of the solvent to the Soret coefficient, measured at ∇T = 0.12, is F s = −0.101 ± 0.003, which is of the order of k B ∇T at infinite dilution [18]. The contribution to the Soret coefficient is 0.343 ± 0.004. Our thermophoretic forces change linearly with the magnitude of the thermal gradient and with the number of beads N (see Fig. 2). This result is consistent with the equation for F T and it shows that all the systems considered here are in the linear regime.
We show in Fig. 3 the simulation data for the Soret coefficient as a function of the type of colloid (symmetric and asymmetric) and its mass. In both cases the Soret coefficient is positive, meaning that the colloid is thermophobic and drifts towards the cold region. In the symmetric case the Soret coefficient increases monotonically with the total mass of the colloid and saturates at high masses. The dependence with mass and the saturation of the Soret coefficient can be modeled using the expression S T = C N (m n − m s )/(m n + m s ) + D N , where C N and D N vary linearly with the colloid mass (see Fig. 3) and m n and m s represent the masses of the colloid and the solvent, respectively. This relationship resembles the kinetic theory equation for mixtures of elastic spheres interacting through r −12 potentials [20], which is the potential we also employed here. The asymmetric mass distribution results in distinctive differences in the Soret coefficient, which does not change monotonically with the colloid mass. This deviation from the symmetric case is particularly clear for longer colloids (N = 11), for which the Soret coefficient features a maximum at small colloidal masses before converging to Soret coefficients that are clearly smaller than those obtained in the symmetric case. We note that this difference in thermophoretic response is driven by mass asymmetry only and we argue that this physical behavior is induced by the coupling of the inhomogeneous We discuss now how the coupling effect leading to the torque is supported theoretically by nonequilibrium thermodynamics [8]. The change in energy of the colloid due to the orientation is defined by where k E is a force constant that determines the change in energy associated to deviations from the preferred orientation of the colloid and n is the unit vector defining its orientation (see Fig. 1). (See Ref. [8] for an example of an orientational vector in a different system, corresponding to a molecular fluid.) The energy above [Eq. (3)] is equal to the rotational kinetic energy of the colloid, where ω is the angular velocity and I is the moment of inertia. The entropy production is defined by where ζ r is the rotational friction coefficient and L nn is a phenomenological coefficient that is defined below. The linear laws describing the fluxes are dn/dt = −L nn k E T −1 n and τ = −ζ r ω, where τ = Iω is the torque. Using Eq. (3) and Eq. (4), we get where I is the moment of inertia. The entropy production associated to the rotation of the colloid in the presence of the thermal gradient is and the corresponding linear flux-force relations, where J q is the heat flux, L αβ ≡ L βα for α = β is the cross phenomenological coefficient, and λ ≡ L qq /T 2 has the usual meaning of the thermal conductivity in the absence of coupling effects, i.e., L nq = L qn = 0. At the stationary state, dn/dt = 0, and the average orientation of the colloid is given by where cos(θ ) = u J q · n, and the orientation is given by the projection of the unit vector, n, along the unit vector defining the direction of the heat flux (thermal gradient), u J q . Equation (9) predicts a linear dependence of the stationary value of cos(θ ) with the thermal gradient, ∇T . We show in Fig. 4 that our simulation results follow the linear dependence with ∇T for colloids of different lengths, N . Further, the strength of the orientation for a given thermal gradient and temperature increases with the colloid length. This dependence can be rationalized in terms of the coefficient L nn = ζ r /I , which can be calculated for the shish-kebab model investigated here (see the Appendix). Indeed, L nn decreases as the length of the colloid increases, contributing to an increase of the orientation [see Eq. (9)]. The definition of L nn implies that the orientation increases with the mass asymmetry, i.e., the mass dipole. This prediction agrees with our computations [see Fig. 4(b)]. We have used our simulation results and Eq. (9) to quantify the dependence of L nq /k E with the mass dipole (asymmetry). The cross phenomenological coefficient increases rapidly with the mass dipole and it saturates at large dipole values (see Fig. 4). The coefficient follows similar trends as L nn increasing in magnitude with the colloid size.
Our results show that the stationary orientation is determined to a large extent by the mass asymmetry (mass dipole, moment of inertia) and geometry (friction coefficient) of the colloid. The orientation influences the thermophoretic force, and this explains the nonmonotonic dependence of the Soret coefficient reported in Fig. 3. We have computed the Soret coefficient at fixed orientations with respect to the direction of the heat flux. The Soret coefficient depends strongly on the angle that the colloid makes with the heat flux vector (see Fig. 5) reaching the lowest or highest value when the colloid is fully parallel or perpendicular to the direction of the thermal gradient. This observation agrees with results reported in Ref. [12]. The results reported in Fig. 5 offer a clear explanation to the dependence of the Soret coefficient reported in Fig. 3. For a given thermal gradient the colloid will move in the direction of the heat flux according to a given thermophoretic force. The force increases monotonically with the mass and saturates to a fixed value in the limit of very large masses. When the colloid features mass asymmetry a coupling of the heat flux is possible and the colloid aligns with the thermal gradient. The alignment is stronger for larger asymmetries [cos(θ ) → 1], and this results in a reduction of the thermophoretic force (see Fig. 5). Consequently, the Soret coefficient will be smaller than the one corresponding to a symmetric colloid with the same mass.
Finally, we have calculated the magnitude of the torque associated to the orientation of the asymmetric colloid. To do this we fixed the colloid at specific orientations in the interval (0, π ) relative to the direction of the thermal gradient. The torque is always positive in the interval (0, π ), indicating that the colloid would rotate with the heavy mass towards the hot region. The magnitude of the torque increases with the number of beads, N. As expected, the torque for the symmetric colloid is zero, consistent with the observation that it does not adopt a preferred orientation. When considering SI units, the magnitude of the maximum torque calculated here (see Fig. 4) τ/ ∼ 3 will scale linearly with the characteristic energy scale of the system under investigation. Using an energy of 1 kJ/mol, corresponding to typical dispersion interactions, we get ∼5 pN nm, to be compared with torques achieved experimentally using metallic nanoparticles, 100 pN nm [14].

IV. CONCLUSIONS AND FINAL REMARKS
The physical principle discussed here might offer new avenues to manipulate colloidal suspensions consisting of colloids with mass distributions that are inhomogeneous. The thermal field will induce a transient torque and therefore a rotation in these colloids, which could be observed in experimental studies. Our work indicates that colloids of similar shape and different internal mass distribution would experience different thermophoretic forces. This effect can lead to the fractionation of colloids of different mass composition in a thermal field, due to their different thermophoretic forces. Experimental studies of colloidal suspensions consisting of colloids with different composition would be very helpful to test experimentally the thermophoretic coupling effect discussed in our work. Our theoretical formulation provides clues on what key variables can be modified to tune the orientational effect. Specifically, we expect that increasing the moment of inertia or reducing the friction will enhance the thermal orientation effect. Further work considering colloidcolloid interactions will be needed to get a full picture on the thermophoretic torque. One key message from our work is that indeed it is possible to induce orientation in axially symmetric colloids using homogeneous fields. This can be achieved by tuning the internal mass distribution of the colloid and, therefore, this physical effect is driven by the internal degrees of freedom of the colloids.

ACKNOWLEDGMENTS
We acknowledge the EPSRC-UK (Grant No. EP/J003859/1) and the EU NanoHeal ITN project grant agreement No. 642976 for financial support. J.M.R. thanks The Leverhulme Trust for the award of a Leverhulme Professorship to visit the Department of Chemistry at Imperial College London. We thank the Imperial College High Performance Computing Service for providing computational resources.

APPENDIX
For the shish-kebab model studied here the friction is given by ζ r = πηL 3 /[3 ln(L/2σ )] [21], where η is the viscosity and L = Nσ is the length of the colloid. For the model used in the paper to describe the mass asymmetry we can define the moment of inertia as where m and m end are the masses of the normal and heavy beads, respectively. r j are the coordinates of the beads along the colloid axis and r COM is the coordinate of the center of mass of the colloid, which is given by The moment of inertia can be simplified to give