Thermomolecular Orientation of Nonpolar Fluids

Frank Römer, Fernando Bresme,* Jordan Muscatello, Dick Bedeaux, and J. Miguel Rubı́ Department of Chemistry, Imperial College London, SW7 2AZ, London United Kingdom Department of Chemistry, Norwegian University of Science and Technology, Trondheim, Norway Process and Energy Laboratory, Delft University of Technology, Leeghwaterstraat 44, 2628 CA Delft, The Netherlands Departamento de Fı́sica Fonamental, Facultat de Fı́sica, Universitat de Barcelona, Martı́ i Franquès 1, E-08028 Barcelona, Spain (Received 25 November 2011; published 7 March 2012)

Thermal gradients are responsible for a wide range of nonequilibrium effects, electron transport (thermoelectricity) [1], mass transport in suspensions (thermophoresis) [2][3][4][5][6], mass separation in liquid mixtures [7][8][9], and nucleation and growth of colloidal crystals [10].Recently it has been shown that temperature gradients can induce orientation in polar fluids, a physical effect that is supported by nonequilibrium thermodynamics theory (NET) and that can be explained in terms of the coupling of a polarization field and a temperature gradient [11].The polarization of the fluid results in sizable electrostatic fields whose strength scales linearly with the temperature gradient.
In this Letter, we show that temperature gradients can also induce the molecular orientation of nonpolar fluids.We call this effect thermomolecular orientation (TMO).The physical origin of this effect cannot be discussed in terms of the polarization-field-heat-flux coupling, although it is important to note that there is not a principle within nonequilibrium thermodynamics theory that precludes observing molecular orientation in nonpolar fluids.
To investigate this hypothesis we use boundary driven nonequilibrium molecular dynamics (NEMD) simulations.NEMD has been successfully employed before to investigate thermodynamic and transport properties of atomic and molecular fluids as well as simple ionic liquids [12][13][14][15].In this work we focus on diatomic fluids.The advantage of using such a model is that it provides simplicity, the necessary ''flexibility'' to change the molecular anisotropy in a controlled way, and as we will see below, it is possible to establish a direct connection with binary mixtures, making it possible to investigate correlations between TMO and the Soret effect [16,17].
We have investigated the TMO effect using a diatomic molecule consisting of two tangent spheres of diameters 1 , 2 , masses m 1 , m 2 , and bond length The bond length was kept constant using a rigid bond through the Rattle algorithm [18].The interaction between the particles is completely repulsive, U ij ðrÞ ¼ 4"ð ij =rÞ 12 , where " is the interaction strength, which we use to define the reduced temperature, T Ã ¼ k B T=", and ij ¼ ð i þ j Þ=2 is defined in terms of the diameters of sites i and j.We note that preliminary simulations indicate that the TMO effect is also present when the intramolecular bond is modeled using a harmonic potential and/or dispersion interactions.The interactions were truncated using a long cutoff of 5, to ensure the discontinuities in the forces between the atoms are small.The tangent sphere model can be seen as a limiting case of a binary mixture, when two atoms of different species are bonded together.We find that the model provides a natural link with the Soret effect in binary mixtures, a connection that we exploit below.The simulations of the binary mixtures were performed using the same interatomic model and cutoff discussed above.In this case the rigid bond was not included.
To investigate the impact of the molecular anisotropy on the thermomolecular orientation we varied the diameter ratio of the two sites in the interval, 2 = 1 ¼ 1, 3=4, 2=3, 1=2 and 1=3 (see Fig. 1).In order to compare the results between the different fluids we fixed the average packing fraction to ¼ ðN=VÞV mol ¼ 0:37, where N is the number of molecules, V is the total volume of the simulation box, and is the molecular volume.The simulation of the binary mixtures involved the same diameter ratios discussed above.
The nonequilibrium molecular dynamics simulations were performed using a modified version of the heat exchange algorithm [12,13] using a bespoke code [14].All simulations were performed using N ¼ 768 molecules in an orthorhombic periodic simulation cell with dimensions ðL x ; L y ; L z Þ ¼ ð3; 1; 1ÞðL=Þ 3 , where L= was varied between 9 and 10.8 to obtain the same packing fraction for each diameter ratio.The heat flux and thermal gradient were generated by scaling the velocities of the molecules in cold and hot layers located at the center and edges of the simulation box, respectively, keeping the total momentum of the box equal to zero.The pressure in the NEMD simulation is constant along the box, ensuring mechanical equilibrium.The equations of motion were integrated using the velocity Verlet algorithm [19] with a reduced time Using parameters for a diatomic molecule such as CO, "=k B ¼ 42:3 K and 0 ¼ 3:2717 A, this corresponds to t ¼ 4 fs.All the simulations involved a preequilibration period of 10 5 time steps at constant temperature, ð2-5Þ Â 10 6 steps to reach the stationary conditions, and a further ð2:5-5Þ Â 10 6 of production.
To quantify the molecular orientation induced by the thermal gradient, we calculated the order parameter, hcos x i ¼ h ũ Á xi, where ũ is the unit vector defining the orientation of the diatomic molecules (see Fig. 1 for the criterion used in this work), x is the vector defining the direction of the gradient in the x direction, and is the angle between both vectors.The order parameter, hcos x i, is zero when the molecular axis does not have a preferred orientation with respect to the gradient along the x axis, and 1 or À1 when the molecular axis is fully aligned parallel/antiparallel to the direction of the gradient.The average temperature in all our simulations was set to hT Ã i ¼ 1:3 or 1.0.The average temperature gradients calculated in the middle region between the hot and cold thermostats, ignoring the thermostatting layer regions, were rT Ã ¼ 0:08, 0.12 and 0.14, which using the parameters given above correspond to, % 1, 1.5 and 1:8 K= A, respectively.
Figure 2(a) shows the order parameter for diatomic molecules with varying size ratio 2 = 1 .We have represented the orientation as a function of the local temperature or local packing fraction in the simulation box.The points represent local properties obtained by dividing the simulation box into 120 layers.We then exploit the symmetry of the box to obtain a final average using both the left and right parts of the simulation box.Our results clearly show that the molecules adopt a preferred orientation as a response to the thermal gradient.The sign of the order parameter indicates that the smaller atom points preferentially towards the cold region.Crucially, the TMO effect becomes stronger with increasing molecular anisotropy [see Fig. 2(a)].The TMO is also a strong variable of temperature and packing fraction.The orientation strongly increases at lower temperatures.This trend is followed by all the molecules investigated here, except the symmetric case, 2 = 1 ¼ 1, where the average orientation is zero.Figure 2(b) shows the orientation in terms of the packing fraction and temperature gradient strength.To isolate these two effects we have focussed on the system 2 = 1 ¼ 1=2.The orientation strongly increases with the packing fraction (lower temperatures).Interestingly, the TMO effect increases with the temperature gradient strength.An analysis of our results for the three different gradients, shows that the orientation scales linearly with the temperature gradient.This is an important result, as it shows the nonequilibrium origin of the TMO.As we will see below this linear dependence can be rationalized using NET.The analysis reported above shows that the TMO response can be tuned using a range of variables; the thermodynamic state (temperature and packing fraction), the molecular anisotropy, 2 = 1 , and the temperature gradient strength.This offers a wide range of possibilities to maximize the effect in practical situations.
We show in the following that the molecular orientation discussed above might be correlated with the Soret or thermophoretic effect, whereby particles of different size and/or mass tend to migrate towards cold or hot sources.To establish this connection we performed simulations of equimolar binary mixtures with size ratios corresponding to those discussed above for the diatomic molecules.The thermodynamic states, average temperature, and packing fraction were the same as the ones investigated above.
The Soret coefficient was computed using [16], s T ¼ ½À1=ðx 1 x 2 Þðrx 1 =rTÞ, where x i is the molar fraction of species i and rx 1 the average molar fraction gradient.We note that in this work s T < 0 indicates that the smaller component migrates towards the cold region.
To clearly show the correlation between the orientation and the Soret coefficient, we have represented the orientation of the molecules versus the Soret coefficient of the corresponding binary mixtures (see in Fig. 3).For a given molecule we represent the average orientation in the simulation box.For m 1 ¼ m 2 and 2 = 1 Þ 1 the Soret coefficient is negative, indicating the smaller site concentrates in the cold region.This observation is compatible with previous work on Lennard-Jones binary mixtures [20].The sign of the Soret effect is compatible with the orientation observed in the diatomic molecules for m 1 ¼ m 2 , i.e., the small site points towards the cold region.
The connection between orientation and the Soret effect is a very important notion as it provides a microscopic principle to rationalize the preferred orientation of the molecules.We can make this connection more explicit by exploring other variables that can drive the Soret effect.One of these variables is the mass [13,20,21].Following the connection between orientation and the Soret effect discussed above, a diatomic fluid consisting of sites with the same size ( should also exhibit TMO.To test this hypothesis we performed additional simulations of both diatomic molecules and binary mixtures for a range of mass ratios, m 2 =m 1 ¼ 1 À 1=20, while the diameter of the particles was 1 ¼ 2 .The results (see quadrant s T > 0 and cos < 0 in Fig. 3) show again a correlation between the Soret coefficient and the molecular orientation.The stronger the separation in the binary mixture the stronger the orientation in the molecular fluid.The sign of the Soret coefficient indicates that the heavier site tends to accumulate at the cold region.This observation correlates well with the net orientation of the molecules, as the heavier site tends to point towards the cold region as well.The main conclusion from our analysis is that the Soret coefficient for binary mixtures provides an approach to predict the preferred orientation of our molecular fluids under thermal gradients.
The thermomolecular orientation discussed above is a new physical effect that has not been considered before, not even in the context of nonequilibrium thermodynamics theory.Hence, in the following we derive the linear fluxforce equations for this new nonequilibrium problem.
We first note that under equilibrium conditions the molecules do not show a preferred orientation.Because of a temperature gradient, as shown in the simulations, the molecules become aligned.We can quantify the degree of alignment through an average orientation vector, h ñi, with ñ ¼ ðcos; sinÞ.To derive the entropy production we consider the energy required to obtain this alignment.Assuming a harmonic potential, E ¼ Ch ñi 2 =2, the force constant C will determine how easy it is to maintain a preferred orientation of the molecules with respect to the isotropic case.Now we derive the entropy production.We consider first the time variation of the energy, @E=@t ¼ Ch ñi Á ð@h ñi=@tÞ.The entropy production is given by, where Jq is the heat flux.The corresponding linear fluxforce equations are given by: where L ij are the phenomenological coefficients, which fulfill the Onsager symmetry relations, L nq ¼ L qn .Solving for the stationary state we get When the gradient is parallel to the x axis, as is the case in our simulations, we get: In NET the ratio L nq =ðCL nn Þ can depend on the temperature but not on the gradient of the temperature.Hence, according to Eq. ( 5) the orientation depends linearly on the temperature gradient.
The inset in Fig. 2(a) shows that the product Thcos x i is constant, within the statistical uncertainty of our computations, for the temperature interval considered here.We therefore find that the ratio L nq =ðCL nn Þ is in good approximation independent of the temperature in this range.Thcos x i increases with the molecular anisotropy, showing the ratio L nq =ðCL nn Þ is determined by the molecular geometry.Because the ratio L nq =ðCL nn Þ is approximately constant for the range of temperatures investigated here, Eq. ( 5) predicts that for a given molecular anisotropy, the orientation curves in Fig. 2 should collapse into a single curve when they are scaled by the corresponding temperature gradient.This notion is confirmed by our results [see inset Fig. 2(b)] showing that the NET captures the correct physical behavior observed in our simulations.
In summary, we have shown that the application of a temperature gradient to a fluid consisting of nonpolar tangent spheres, which is representative of a molecular fluid consisting of diatomic molecules, induces the alignment of the molecular axis along the temperature gradient.This TMO effect scales with the molecular asymmetry, being stronger for molecules featuring sites with different diameter and/or mass.By inspecting the dependence of the orientational response with the size and mass of the molecular sites we have been able to establish a correlation between the TMO effect and the Soret effect.This is an important physical principle that we believe can help in the future to control the strength of the effect by tuning the molecular architecture.The strength of the orientation observed using the present model may appear small.However, it is of the same order as the orientation observed in polar fluids under temperature gradients, which results in significant electrostatic fields, $10 4 -10 6 V=m for thermal gradients, 10 6 -10 8 K=m, which are achievable in experimental setups [5,22].Interestingly, the orientation is also similar in magnitude to that observed in water at the liquid-water interface [23].The TMO effect should be observed at all length scales provided the gradient is strong enough.The NET supports this notion.Since the model investigated in this work is also applicable to colloidal dumbbells, we expect that experiments at the m scale using optical techniques may provide a route to directly observe the physical effect discussed here.
We would like to thank the Leverhulme Trust for financial support and the Imperial College High Performance Computing Service for providing computational resources.F. B. would like to thank the EPSRC for support in the form of a Leadership Fellowship (EP/J003859/1).*f.bresme@imperial.ac.uk

FIG. 2 (
FIG. 2 (color online).(a) Orientation of the diatomic moleculeversus the reduced temperature for different size ratios.From bottom to top 2 = 1 ¼ 1, 3=4, and 1=3.The mass in all cases is m 2 =m 1 ¼ 1 and the temperature gradient rT Ã ¼ 0:08.The inset shows the same plot but for the product hcosð x ÞiT Ã in the y axis.This plot is used to test the NET equations derived in the Letter.See main text for details and discussion.(b) Dependence of the orientation with particle packing fraction for 2 = 1 ¼ 1=2 as a function of the temperature gradient, rT Ã ¼ 0:08, 0.12 and 0.14.The average temperature for all the systems is hT Ã i ¼ 1:3.The inset shows the rescaled plot for the orientation, i.e., hcosð x Þi=ðrT Ã =0:08Þ, to test the NET equations derived in the Letter.See discussion in the main text for details.