A dissipative particle dynamics study of the realignment of a nanodroplet of a nematic in a weak external magnetic field

Abstract. We present a dissipative particle dynamics (DPD) approach for simulating the realignment of a nematic nanodroplet suspended in an isotropic fluid following a switch in the direction of an applied external magnetic field. The interaction of the mesogens with the external field is weak relative to the inter-molecular interactions. The simulations were used to investigate the way orientational equilibrium is re-established. The results reveal that the realignment process of the nanodroplet is consistent with its fluid structure. The reorientation of the nanodroplet as a whole is found to be caused by an internal structural rearrangement rather than a coherent rotation of the centres of mass of the mesogens about the centre of the nanodroplet. The switch in the field direction furthermore is found to induce a transient spatial variation in the orientational order of the long axes of the mesogens: the orientational order parameters decreases on moving from the core of the nanodroplet to the surface in contact with the isotropic environment. The results highlight differences between the time evolution of the orientation of the long molecular axes in the field and the rotations of the centres of mass of the mesogens about the centre of the nanodroplet.


Introduction
Dissipative particle dynamics (DPD) [1][2][3][4][5] is a computational method based on the concept of "matter particles", each representing a cluster of atoms or mesogens, interacting via soft potentials and subject to dissipative and fluctuating forces. Thus polymeric systems can be described by simple particle-spring models with long chains represented by 20 to 50 particles [4]. The technique was originally developed for describing isotropic phases with the important advantage of providing a way of describing the qualitative behaviour of a fluid system at relatively long times. A dissipative particle dynamics description does not provide a picture of the static and dynamic properties of the oriented fluid phases in terms of molecular properties. However, it affords a rationalisation of complex dynamic processes in ordered liquid phases in physical terms, which can be complementary to more detailed, atomic-based simulation tools. This is expected to be particularly true when considering classical rheological experiments monitoring the distribution of the director orientation patterns a e-mail: y.k.levine@phys.uu.nl (corresponding author) b e-mail: antonino.polimeno@unipd.it in nematic and smectic phases in the presence of mechanical or field-induced stresses, which are known to exhibit dynamic phenomena on different time scales, and moreover to be highly dependent on the system geometry, order parameter as well as field intensity.
We have demonstrated [6] that orientational and translational order can be induced with the "soft" interaction potential employed in conventional DPD algorithms. The systems used consisted of chains of particles whose flexibility has been hindered by an externally imposed harmonic potential. The computational effort is relatively modest compared to that demanded by other simulation techniques, thus establishing DPD as a convenient coarsegraining technique for numerical simulations of liquidcrystalline phases.
We have here used the DPD approach to simulate a non-equilibrium situation, namely the realignment of a nanodroplet of a nematic liquid crystal suspended in an isotropic fluid under the effect of an external magnetic field. Experiments utilizing such systems have been reported in [7][8][9]. The interaction of the mesogens with the external field is weak relative to the inter-molecular interactions.
The paper deals separately with the rotations of the centres of mass of the mesogens about the centre of the nanodroplet and realignment of the long axes about the centres of the mesogens. We furthermore investigate the spatial inhomogeneities in the properties of the nanodroplet during the realignment process. We find that the overall orientation of the nanodroplet in the field is consistent with a fluid structure and caused by an internal structural rearrangement rather than a coherent rotation of the centres of mass of the mesogens about its centre. The switch in the field direction furthermore induces a transient spatial variation in the orientational order of the long axes of the mesogens: the orientational-order parameters decreases on moving from the core of the nanodroplet to the surface in contact with the isotropic environment. Nonetheless, realignement process of the sample-averaged orientational-order parameters can be well described by the behaviour of a reorienting rigid monodomain nematic sample.
2 Dissipative particle dynamics simulations of ordered phases 2.1 Dissipative particle dynamics of nematic mesogens The DPD methodology for the simulation of orientationally ordered phases has been described in detail previously [1][2][3][4][5] and will only be summarised here. The equations of motions governing the behaviour of particles of identical mass are where r i and v i are the position and velocity of the i-th particle, respectively. The conservative force, f C , arises from a pair-wise potential, while f D and f R represent the effective dissipative and fluctuating forces, respectively. The potential giving rise to the conservative forces acting on identical particles has the form of a sticky soft interaction: where a is the force parameter, setting the overall repulsion strength between particles, ε rep ij and ε att ij modulate the repulsive and attractive interactions, respectively, r c is the cut-off radius of the repulsive interaction and r ij the distance between the particles. The attractive interaction corresponding to a sticky potential operates in a thin shell of width δr around the particle. We have here set δr = 0.1r c . Our choice for the form of the attractive interaction has been made in the spirit of the DPD framework Reduced energy though other forms can be used [4]. The conservative force f C acting on the i-th particle is then simply given by wherer ij is the unit inter-particle. The dissipative force f D is defined as a quantity proportional to the relative velocity between the particles, Here v ij is the projection of the relative velocity on the vectorr ij and w C (r ij ) = 1 − r ij /r c when r ij < r c + δr and zero elsewhere. The dissipation strength is determined by γ, and the corresponding effective random force is defined as f R = σξ(t)w C (r ij )r ij / √ dt, where ξ(t) is a white-noise random term. A fluctuation-dissipation relation holds for the parameters of the dissipative and random forces, In the simulations below we make use of reduced units where r c , m and a are the units of length, mass and force, respectively. The other units are defined in Table 1.
Each mesogenic unit in the simulation is a chain of n particles linked by elastic forces given by the classical spring force f i = −k * (r ij − r eq )r ij = −f j , with k * the reduced spring constant, r ij the distance between the two joined particles and r eq the equilibrium particle-particle distance. In addition, an elastic force with an elastic constant k * e is introduced between the first and last particles: This force primarily determines the linear extent of the mesogenic unit and random coils are obtained when it is omitted.
A spherical simulation box of radius R centred at the origin was used in this study. The mesogenic units were confined to the box by the application of a shortrange external radial force to each particle, i: f i,sphere = k * in exp[−2.5(R−r i )]r i , wherer i is the unit position vector of the particle with respect to the origin of the simulation sphere.
The equations of motion are solved numerically using the standard Verlet algorithm coupled with a thermostat rescaling the velocities. Periodic boundary conditions were not applied as we are here dealing with a single nanodroplet in the simulation sphere. A linked-cell algorithm [10] was implemented.

Realignment of nematics in an external field
We shall now consider the alignment of the model mesogens in the nandroplet under the influence of an external field. The applied field acts on the two terminal particles of the mesogenic unit (i.e. on the end-to-end vector of the mesogen). The interaction used here corresponds to the application of an external magnetic field and takes the form of an interaction between a field and an induced dipole: where ε * ext is the strength of the interaction and ν i (t * ) is the unit end-to-end vector of the i-th mesogen and H ext is the unit vector of the external field. In the work presented here, the external field lies in the XZ plane of the simulation sphere. The calculation of the forces from the potential is described in the appendix.

Characterisation of the nanodroplet shape
The shape of the nanodroplet is readily characterised by evaluating its tensor of inertia using the Cartesian components of the position vector of the centres of mass, r cm i , of the N mesogens relative to the centre of mass of the nanodroplet. The principal axes of inertia and their orientation in the laboratory are obtained from a diagonalisation of the matrix I αβ .

Orientation of the long axes of the mesogens in the simulation sphere
The rotational anisotropy of the long axes of the mesogens with respect to the fixed Z-axis is characterised by the orientational-order parameter S zZ (t * ): where θ i (t * ) is the angle made by the end-to-end vector of the mesogen with the Z-axis at time t * . Similarly, we calculate the order parameters S zX (t * ) and S zY (t * ) of the end-to-end-vector relative to the X-and Y -axes, respectively. As the three parameters are not independent, but related by S zX (t * ) + S zY (t * ) + S zZ (t * ) = 0, we shall here only discuss the behaviour of the latter two order parameters.

Rotation of the centres of mass of the mesogens about the centre of the nandroplet during realignment
The whole-body rotation of the mesogens as the result of the switching of the applied field can be followed by calculating the angular excursion of their centres of mass about the three axes of the simulation sphere, θ X , θ Y and θ Z . This can be achieved in two ways. The angles can be computed from the spherical components of the position vectors relative to the axes of the simulation sphere. This, however, requires numerical checks to take account of sudden changes in the angle of azimuth across the 0 • /360 • boundary. An equivalent, though simpler procedure which does not suffer from this makes use of the projections of the position vectors of the centre of mass on to the ZY , XZ and XY planes, respectively. The rotations are calculated by evaluating the vector products between the components at times t * and t * +dt * and making use of the fact that we are dealing with small angles so that sin θ ≈ θ.
In this way we obtain the following expressions for time evolution of the rotation angles about the three fixed axes, X, Y and Z during the course of the simulation: The angled brackets denote an average over the mesogens in the nanodroplet.

Results of the simulations 3.1 System parametrisation
The internal conformation of a model mesogenic unit is characterized by 4 parameters: the number of particles n, the two elastic spring constants, k * and k * e and the equilibrium distance between joined particles, r eq . We have previously shown that on using identical mesogenic units each containing 10 particles with r eq = r c , ε * rep ij = 24, ε * att ij = 0, k * = 100 and k * e = 500, at a reduced density ρ * = 0.8 the system forms a nematic phase in the reduced temperature range 0.25 ≤ T * ≤ 0.5 [6].
The current simulations make use of a nanodroplet of a nematic phase suspended in an isotropic environment. The isotropic environment consists of mesogenic units in a random-coil configuration, that is for these mesogens k * e = 0. We found that this prescription avoided the numerical instabilities reported in [5] for simple fluids of spheres. The radius of the nanodroplet was taken to be 0.7R, with R the radius of the spherical simulation box. On taking R = 143, the system configuration consisted of 36786 stretched mesogens in the nanodroplet and 69684 random coils forming the isotropic environment. Under this condition the external confining force, k * in = 24, only acts on the random coils at every stage of the simulation run.
Time steps dt * = 0.04 were used and the integration scheme yielded temperature fluctuations of less than 0.25% in the equilibrated state. A thermostat was applied every 5 time steps during the simulation run.
We were unable to produce a stable nanodroplet using purely repulsive potentials and therefore described the particles as sticky spheres. To this end, the strengths of the inter-particle interactions between mesogens of the same configuration were taken to be ε * rep ij = 24 and ε * att ij = 60, while those between particles in the stretched and randomcoil configurations (in the nematic and isotropic phases, respectively) were ε * rep ij = 60 and ε * att ij = 0. This choice of parameters produced sufficient surface tension between the nanodroplet and its isotropic environment to suppress the exchange of mesogens between the two phases. We emphasise, however, that we did not carry out an exhaustive search of parameter space and that consequently other combinations of values may prove to be as effective.
The system configuration produced in this way was found to be stable in the reduced temperature range 0.25 ≤ T * ≤ 0.45. The nanodroplet underwent a nematicisotropic transition at the upper limit of this interval. A typical image of the system at T * = 0.325 is shown in Figure 1. The image sets out the spatial positions of the terminal particles of 15000 mesogens selected at random. The particles are connected by a virtual bond and the image shows a view projected on the XZ plane of the simulation sphere.
Realignment experiments on nematic phases are usually carried out with magnetic field of up to about 10 T [11,12]. Under these conditions the interaction of the nematogens with the external field is much weaker than the inter-molecular interactions in the system. We shall here emulate this situation and choose ε * ext to be as small as practically possible. We have found that on taking ε * ext = 0.1, the orientational energy in the applied field amounts to about 3% of the inter-molecular energy. With this value of ε * ext we find that the equilibrium value of S zZ increases from 0.76 to 0.84 on applying the orienting field, indicating that the external interaction only causes a small change in the intrinsic system properties. This change was found to be due entirely to an alignment of the end-to-end vectors of the mesogens in the field. The dimensions of the mesogens are not affected by the presence of the external field. We note here that the internal elastic energy of a model mesogen is 3 orders of magnitudes larger than the orientational interaction. In order to check whether the realignment processes depend on the value of the field, we carried out a single simulation run with ε * ext = 0.02, 5 times smaller than that used in the results presented below. Now the equilibrium value of S zZ increased from 0.76 to 0.775 on applying the orienting field.
The starting point of the simulations was a thermally equilibrated system configuration of mesogens with the external field applied along the Z-axis of the simulation sphere. We note that the field only acts on the stretched mesogens constituting the nanodroplet.
At time t * = 0, the external field was switched in the XZ plane to a different orientation at an angle θ ext relative to the Z-axis. The properties of the system were then monitored for t * ≤ 5000 during the return to orientational equilibrium along the direction of the field. The The results of the simulations reported below for a switch angle of 60 • were virtually identical with those obtained with ε * ext = 0.02 provided the time axis was scaled by a factor of 5 as expected theoretically. This is a strong indication that the mechanism of alignment is independent of the field value for our choice of parameters.
Finally, it is important to emphasise that the latter case, a switch of the applied magnetic field by 90 o , is a special case [12]. The field-mesogen interaction, equation (5), has an extremum when the angle between the field and the long axis of a mesogen is 90 o . Consequently, the external torque working on the mesogen works in a clockwise direction when the angle is less than 90 o and counter-clockwise when it is greater than 90 o . The torque changes direction as the long axis of the mesogen swings through 90 o as the result of orientational fluctuations and thus enhances the effect of the fluctuations. This degeneracy should produce two counter-rotating mesogen populations, resulting in a symmetric distribution of rotations about the Y -axis of the simulation sphere. Furthermore, anomalous changes in the internal structure of the nanodroplet are expected to be observed during the realignment process.

Equilibrium properties
The eigenvalues (principal values) of the moment-ofinertia tensor of the nanodroplet were calculated to be 386.5, 393.4 and 460.2 about the x-, y-and z-axis of the nanodroplet, Figures 2. The equilibrium shape of the nanodroplet is thus essentially a prolate sphere with its long axis oriented along the Z-axis of the simulation sphere. The spherical nanodroplet was found to be spatially homogeneous except for an outer shell in contact with the isotropic environment. The heterogeneity in the properties of the nanodroplet can be monitored by dividing it into a series of shells. To this end we calculated the position vector of the centre of mass of each mesogen, The corresponding values of the order parameter S zZ were found to be 0.865, 0.871, 0.867 and 0.780. The latter value obtained for the outer shell is lower than the rest as a result of the constraints on the alignment of the long axes of the mesogens by the curved surface of the nanodroplet.

Shape of the nanodroplet during alignment
The principal values of the moment-of-inertia tensor of the nanodroplet vary smoothly by around 15% during the realignment process following the switching of the magnetic field by less than 90 o , Figures 2. However, anomalous changes are observed following a 90 o switch in the direction of the field, as expected from the degeneracy of the field-mesogen interactions immediately after the switching of the field. Importantly, we note that the nanodroplet does not reach an equilibrium shape during the simulation run.
It can be seen that the limiting values of the principal moments of inertia about the x-and z-axes of the nanodroplet reached at long time, t * > 3000, are significantly different from their initial, equilibrium, values. This indi-cates that the internal configuration of the nanodroplet changes irreversibly during the realignment process. Nevertheless, the principal moment of inertia about the y-axis of the drop is essentially constant, and this is consistent with the switching of the field in the ZX plane of the simulation sphere. Importantly, we find that the magnitude of the change increases on increasing the switching angle of the field, Figures 2. It is important to note that in our simulations the time scale of translational diffusion is comparable to that of the rotational realignment processes. We believe that the change of the nanodroplet shape during the realignment process can be ascribed to the combination of the random distribution of the centre of mass of the mesogens and their anisotropic translational diffusion.
An estimate of the magnitude of the whole-body rotation of the nanodroplet in the XZ plane away from the Z-axis of the simulation sphere towards the field may be obtained from the time dependence of the eigenvector corresponding to the long axis of the nanodroplet, Figure 2c, on assuming that the changes of shape are negligible. The angles of rotation of this eigenvector about the Y -axis of the simulation sphere obtained in this way are shown in Figure 3 for all the cases studied. The limiting values reached at long times are found to be close to the switching angle applied.
We shall now address the question as to whether the nanodroplet rotates as a rigid entity and if the changes of the overall shape found here are due to spatial inhomogeneities in its configuration. The averaged rotations of the centres of mass of the mesogens in the core of the nanodroplet and the three spherical shells about the Y -axis of the simulation sphere are shown in Figure 5. It can be clearly seen that the nanodroplet exhibits a gradient of rotational motions: mesogens at the core undergo significantly smaller rotations about the Y -axis than mesogens in the outer shell. Moreover, the average rotation of the centres of mass of the mesogens about the centre of the nanodroplet is considerably smaller than the switching angle of the applied field in every region. This strongly indicates that the rotation of the long axis of the moment-of-inertia tensor of the nanodroplet, Figure 3, occurs as the result of a structural rearrangement of mesogens in the nanodroplet rather than a coherent rotation of the centre of mass of the mesogens characteristic of the rotation of a rigid nanodroplet.

Orientational distributions of the long axes of the mesogens during realignment
The normalised orientational probability distribution functions of the long molecular axes of the mesogens immediately following the switch in the direction of the applied field and those found at the end of the simulation time t * = 5000 are shown in Figures 6a and b, respectively. The initial distribution functions are clearly broad and unimodal. However, only the distribution functions for the 75 o and 90 o switching angles exhibit a significant probability of orientations of the long axes of the mesogens at angles larger than 90 o relative to the field. The final distribution functions are also unimodal, Figure 6b, and identical for all cases except that of a 90 o switch of an applied field. It thus appears that the width of the initial orientational distribution functions and the fluctuations of the orientations of long axes of the mesogens induced by the DPD algorithm used here are not sufficient to tilt a significant population of molecular axes to angles larger than 90 o for switching angles θ ext ≤ 75 o . Consequently, two counter-rotating orientational populations are only observed for a switching angle of 90 o . Moreover, we note that the normalised orientational probability distribution functions shown in Figure 6b are virtually identical with those observed at equilibrium.
The realignment process of the long axes of the mesogens in the nanodroplet can be monitored in more detail by calculating the temporal behaviour of the orientation order parameters S zY (t * ) and S zZ (t * ) following a switch in the orientation of the applied field. The order parameter S zY (t * ) monitors the excursions of the long axes of the mesogens out of the XZ plane during the realignment process. The order parameter S zZ (t * ) on the other hand, monitors the reorientation of the distribution of mesogen axes from the Z-axis at time t * = 0 towards the field. Again, we shall evaluate the order parameters separately for the core of the nanodroplet and the three outer shells.
The behaviour of the order parameters S zY (t * ) are shown in Figure 7. It can be seen that for all the switching  Figure 8. It can be seen that the switching of the applied field induces a significant spatial gradient in the orientational order from the core to the outer shell during the realignment process. Nonetheless, a homogeneous distribution is recovered at long times, t * > 4000 except for the cases of switches of the applied field by 30 o and 90 o , where the inhomogeneity persists at the end of the simulation run. For all other angles the values of the order parameters of the long axes of the mesogens calculated at the end of the simulation run with respect to the field are close to the initial equilibrium values. Thus, an orientational equilibrium is re-established along the applied field. This is in contrast to the mechanical properties described in Section 3.4.
The spatial inhomogeneity in the orientational order found in our simulations is not revealed experimentally, where commonly order parameters averaged over all the mesogens in the sample are monitored. We shall therefore now consider the behaviour of the sample-averaged order parameters and address the question as to whether they can be described by the analytical approaches used to describe experimental data. For simplicity we shall only consider the torque-balance model [13,14] for the realignment of a rigid monodomain nematic sample in a magnetic field. The equation of motion of the director of the sample is given by where β(t) is the angle made between the director and the external field at time t following an initial switch by an angle β 0 at time t = 0. We shall now assume that the orientational distribution of the long axes of the mesogens retains its initial unimodal form during the realignment process. Under this assumption we can express the orientation of any mesogen in the nanodroplet in terms of two successive rotational transformations: i) a rotation from the fixed Z-axis to the instantaneous axis of symmetry of the distribution and ii) a rotation from the instantaneous symmetry axis to the long axis of the mesogen. Transformation i) can then be used to characterise the behaviour of the simulated system as a monodomain. The transformation properties of orientational-order parameters under rotations have been discussed in detail in the literature [15,16] and it is known that the order parameter S zZ (t * ) calculated in the laboratory frame for a uniaxial distribution can be written as where S zZ (t * ) 0 is the intrinsic order parameter of the distribution calculated relative to the symmetry axis, and θ A (t * ), is the angle made by the symmetry axis with the fixed Z-axis at time t * . On setting S zZ (t * ) 0 = S zZ (t * = 0), we can extract the order parameter of the instantaneous orientation of the symmetry axis of the distribution from the normalised order parameter S z (t * )/S z (t * = 0). The angle β(t * ), equation (9), between the director of the nematic monodomain and the applied magnetic field can be calculated since S zZ (t * )/S zZ (t * = 0) = 0.5 3 cos 2 [β 0 − β(t * )] − 1 . (11) The sample-averaged order parameter S zZ (t * ) and the fit of the extracted angle β(t * ) to equation (9) are shown in Figures 9a and b, respectively. Surprisingly, the angles obtained from the simulations exhibit a monoexponential behaviour over large portions of their decays. The decay constants however, vary almost linearly from 840 at θ ext = 30 o to 1095 at θ ext = 75 o . It thus appears that the results presented here for the sample-averaged orientation order of the long molecular axes of the mesogens is described well by the torque-balance equation, even though the detailed behaviour is difficult to reconcile with its underlying assumptions.

Conclusions
We have here presented a dissipative particle dynamics study of the realignment kinetics of a nematic nanodroplet suspended in an isotropic fluid following a switch in the direction of an applied field. The results reveal significant spatial inhomogeneities in the properties of the nanodroplet, consistent with its fluid structure. The reorientation of the nanodroplet as a whole is found to be caused by an internal structural rearrangement rather than a coherent rotation of the centres of mass of the mesogens about the centre of the nanodroplet. The switch in the field direction furthermore is found to induce a transient spatial variation in the orientational order of the long axes of the mesogens: the orientational-order parameters decrease on moving from the core of the nanodroplet to the surface in contact with the isotropic environment. These results highlight the differences in the reorientation of the long molecular axes about the centres of mass of the mesogens, and the rotations of the centres of mass about the centre of the nanodroplet.
This work was financed by the Italian Ministry for Universities and Scientific and Technological Research, project PRIN ex-40% and FIRB RBNE01P4JF.

Appendix A.
We shall here consider the interaction, equation (5), between the unit end-to-end vector of a mesogen, ν = {∆x/r, ∆y/r, ∆z/r} and a magnetic field applied in the XZ plane of the simulation box, H ext = {H x , 0, H z }.
The forces F and −F are now applied to particles 1 and n of the mesogen.