Effect of dipolar interactions on the phase behavior of the Gay–Berne liquid crystal model

Effect of dipolar interactions on the phase behavior of the Gay–Berne liquid crystal model Mohammed Houssa and Luis F. Rull Departamento de Fı ´sica Ato ´mica, Molecular y Nuclear, Area de Fı ´sica Teo ´rica, Universidad de Sevilla, Aptdo 1065, Sevilla 41080, Spain Simon C. McGrother School of Chemical Engineering, North Carolina State University, Raleigh, North Carolina 27695-7905 ~Received 29 May 1998; accepted 21 August 1998! A computer simulation study of the phase behavior of the dipolar Gay–Berne liquid crystal model is presented. The phase transitions are determined with isothermal–isobaric ~NPT!Monte Carlo simulations, utilizing the reaction field method. The electrostatic forces are found to have a considerable effect on the nature of the observed phases, but the density at which the isotropic fluid becomes unstable with respect to partially ordered phases is seen to be remarkably insensitive to the strength of the dipole. We pay particular attention to the structure of the mesophases, combining information from several singlet and pair distribution functions to build up an accurate picture of the molecular arrangement of the systems. © 1998 American Institute of Physics. @S0021-9606~98!51544-9# I. INTRODUCTION It is becoming increasingly evident that the liquid crystalline phase behavior of complex molecules is dictated, principally, by the nonsphericity of the inflexible regions of the constituent molecules. Since the 1940’s, it has been recognized that the simplest liquid crystal ~the nematic!could result from excluded volume effects only. In more recent years, other, more complex, liquid crystalline phases have been reported for model molecules interacting only via anisotropic repulsive potentials. Both smectic1and columnar2 phases have been observed for these simple models. Even the spontaneous polarity of the tilted chiral smectic-C*phase has been rationalized purely in terms of packing effects of the repulsive cores.3 Of course, it is only in these idealized theoretical models that the attractive and repulsive elements of the intermolecular interaction can be considered in isolation. In any real system, both forces will be present. There are, however, colloidal systems ~e.g., tobacco mosaic virus!, in which the nonspherical molecules have a negligible attraction for each other and the phase transitions are thought to be almost entirely driven by excluded volume effects. True to theoretical prediction, such systems do indeed exhibit a number of mesophases at appropriate concentrations and temperatures. But these systems are exceptional; in all thermotropic mesogens, attractive forces are not only present, but frequently comprise complicated components such as multipolar forces and p -electron interactions, in addition to the usual London dispersion effects. It is important, therefore, that the role of the attractive forces in stabilizing the mesophases is investigated. Frenkel and co-workers have shown that hard prolate and oblate ellipsoids of revolution4and hard spherocylinders5have a rich phase diagram with nematic, and in some cases smectic and columnar phases, depending on the molecular size and shape. These findings reinforce the ideas of Onsager,6,7 that the anisotropic shape of the molecules is the principal driving force behind mesophase formation. The attractive forces ~van der Waals, multipolar, dispersion, etc.! are of secondary importance. Despite their lesser role, an investigation of the influence of such attractive forces is necessary if a full understanding of liquid crystalline phase behavior is to be achieved. Probably the simplest idealized liquid crystal ~LC!model which incorporates both repulsive and attractive terms, and consequently the most widely studied, is the Gay–Berne ~GB!potential. Depending upon the choice of parameters this model may be either oblate ~disclike!or prolate ~rodlike!. The majority of real mesogens are rod shaped, thus we focus in this study on the prolate Gay– Berne model. Both the shape- and energetic anisotropy are adjustable. A minimum elongation is required before orientationally ordered fluid phases become thermodynamically stable. Similarly, there is a minimum energetic-anisotropy, below which no spatially inhomogeneous fluid phases are observed. These two minimum values are partially coupled. Obviously there is a vast parameter space to explore; it is necessary, therefore, to limit our investigation to a single energetic anisotropy k 855~ k 8is the ratio of the potential well depths of side-by-side and end-to-end configurations!,a value which yields a rich phase diagram for elongations of interest. Studies of the Gay–Berne potential with k 53~ k is the molecular elongation, and is defined as the ratio of the two principal axes of the ellipsoidal core!and k 855 have proven that the isotropic ~I!, nematic ~N!, smectic-B ~SmB!, and crystalline phases can be formed along different isotherms.8,9 Upon increasing the elongation of the GB particle ~i.e., k .3!, the smectic-A ~SmA!phase appears between the nematic and smectic-B phases, for certain temperatures. For hard ellipsoids with an elongation a/b53~with aand bthe JOURNAL OF CHEMICAL PHYSICS VOLUME 109, NUMBER 21 1 DECEMBER 1998 95290021-9606/98/109(21)/9529/14/$15.00 © 1998 American Institute of Physics Downloaded 24 Sep 2003 to 150.214.138.210. Redistribution subject to AIP license or copyright, see http://ojps.aip.org/jcpo/jcpcr.jsp lengths of the two principal molecular axes!, the isotropic and crystalline phases are separated only by a nematic region. Thus we can confidently state that the existence of the smectic phases is directly related to the attractive forces in the GB model ~the softness of the repulsion is unlikely to cause a qualitative difference in phase behavior!. Studies of dipolar liquid crystal models have also been reported in recent years. Weis, Levesque, and Zarragoichoechea, in a series of papers,10–13 studied how the presence of point dipoles influenced the structure of the LC phases of a variety of hard-core models. Selected state points, representative of various mesophases, were simulated with Monte Carlo, employing the computationally expensive Ewald summation method, principally in the canonical ensemble. Two major conclusions were drawn: central dipoles promote layered structures and terminal dipoles have little effect on the structure of the phase. Little could be concluded as to the stability of the various phases, since only a very limited number of state points were investigated. McGrother et al. constructed full isotherms for several spherocylinder plus dipole models. For central dipoles,14 short-ranged order in the isotropic phase made the I–LC phase transition difficult to pinpoint; indeed, for sufficiently strong dipoles, compression led only to metastable glassy phases and no true LC phases were observed. A much greater degree of hysteresis is noted at the I–N transition than in the hard spherocylinder system. The presence of the central dipole ~both longitudinal14 and transverse15!unambiguously stabilizes the smectic-A phase. This layered phase is seen to be stable at pressures and densities corresponding to less ordered phases in the nonpolar system. When the dipole resides in the terminal position,16 different behavior is noted. Again hysteresis is seen at the I–N transition, but now the smectic-A phase is destabilized, and the N–SmA transition is postponed to higher density and pressure as a result of the electrostatic interaction. In all cases no ferroelectric phases were noted. For the Gay–Berne model, some numerical results for dipolar systems have been offered. Principally Satoh et al.17,18 have used the reaction field method to create several isochores for such systems. For GB molecules ~ k 53 and k 855!with central longitudinal dipoles,17 the I–N transition is again seen to be insensitive to the dipole strength. On the other hand, the N–SmA transition temperature is seen to increase with the dipole moment, i.e., the smectic-A phase is stabilized. The temperature at which crystallization occurs does not vary much with the strength of the multipole. In terms of structure, the layering in the smectic-A phase is observed to be much sharper in the dipolar system, at low temperatures; near the N–SmA transition, the effect is less pronounced. When the longitudinal dipole is placed in the terminal position,18 the I–N transition temperature TIN increases with increasing dipole moment m *5( m 2/«0 s 0 3)1/2 where s 0is the contact distance, and «0is the energy ~ignoring the dipole!of a pair of GB molecules in the side-by- side arrangement. TNS exhibits similar m dependence. The structure of the resultant smectic phase is seen to be bilayered, with a notable degree of interdigitation. Interesting phases were found by Berardi et al.19 when the longitudinal dipole is in either the central or terminal position. This study focused primarily on the smectic phases, and as such large system sizes were employed ~N51000 and 8000!. Three temperatures were selected at a fixed density r *5 rs 0 350.3, corresponding, in the absence of dipolar interactions, to isotropic, nematic, and smectic-A phases, respectively. The addition of point dipoles with m *52 is seen to have no influence on the overall nature of the observed phase, at least at the three chosen temperatures. In both LC phases, the presence of the dipole in either position is seen to increase the orientational order slightly. For central dipoles, the smectic has the usual monolayer structure, but a characteristic split of the second peak in the radial pair distribution function leads Berardi et al. to describe this phase as smectic-B. Shifting the dipole toward one end of the molecule eliminates this hexagonal order within the plane, but several new features are noted. The projection of the pair distribution function shows a splitting in the first peak, which is indicative of a bilayer structure. The first peak is at a distance which suggests significant interdigitation. Particles in each layer are locally ferroelectric ~i.e., the dipoles all point in the same direction!, and the direction of polarization alternates between successive layers. These findings suggest that the phase is a bilayered antiferroelectric smectic phase. However, simulation of huge systems (N58000), and visualization of the phase, shows that the situation is somewhat more complex. In each of the smectic layers, there is a noticeable dislocation, at which the center of the layer shifts by around a molecular length ~due to the toroidal boundaries, there is second dislocation in each layer, which cancels the effect of the first!. Even starting in a perfect smectic-A phase, the system was seen to relax into the eloquently named modulated antiferroelectric bilayer stripe domains. Gwo ´z ´dz ´et al.20 used molecular dynamics ~MD!to study the influence of transverse central dipoles on the phase diagram. These authors take a slightly different approach: rather than cooling along an isochore, they compress along an isotherm, for only one value of the dipole moment ( m * 50.5). They find that the I–SmA transition occurs at the same density with or without the dipolar forces present. The only distinctions between the polar and nonpolar case are the higher degree of pretransitional order in the dipolar system, and better defined positional order in the smectic phase. The authors report smectic-C order in the dipolar case, but the tilt angle is small and the tilt disappears as the density is increased still further. Most recently Houssa et al.9studied a GB model decorated with a central longitudinal dipole. These authors report the complete suppression of the nematic phase for sufficiently strong dipoles. The phase sequence for the nonpolar system is I–N–SmB, but with the inclusion of a dipole of reduced moment m *52.5, the isotropic phase spontaneously orders directly to the smectic-B phase upon compression. Thus from earlier studies of dipolar LC models we can conclude that the dipole causes only a mild perturbation to the initial I–LC transition, which appears to be almost entirely dependent upon the short–range repulsive forces acting between anisotropic molecules. Depending on the position and orientation of the dipole in the molecular frame, the nematic or smectic phases can be promoted at the expense of 9530 J. Chem. Phys., Vol. 109, No. 21, 1 December 1998 Houssa, Rull, and McGrother Downloaded 24 Sep 2003 to 150.214.138.210. Redistribution subject to AIP license or copyright, see http://ojps.aip.org/jcpo/jcpcr.jsp the other. The most evident feature of all previous studies is that the fine detail of the smectic structure is sensitive to the dipole. Interdigitation, bilayers, striped domains, and even tilted phases are reported for various systems. In the next section we describe the potential model, and in Sec. III we give details of the simulation methodologies that have been employed. The results of the simulation studies are presented and discussed in Sec. IV, and we draw our conclusions in the final section. II. POTENTIAL MODEL In our simulations, the interaction energy of two molecules iand jis given by Uij5UGB~r ˆij,u ˆi,u ˆj!1U mm ~r ˆij,u ˆi,u ˆj!,~1! where UGB(r ˆij,u ˆi,u ˆj) is the GB potential21 UGB~r ˆij,u ˆi,u ˆj!54«~r ˆij,u ˆi,u ˆj! 3 F S s 0 rij2 s ~r ˆij,u ˆi,u ˆj!1 s 0 D 12 2 S s 0 rij2 s ~r ˆij,u ˆi,u ˆj!1 s 0 D 6 G ,~2! where u ˆiis the axial vector of molecule iand rij is the distance between the centers of mass of molecules iand j, r ˆij5rij/rij is a unit vector along the intermolecular vector rij5 u ri2rj u , where riand rjare the positions of the centers of mass of molecules i and j, respectively. Here s (r ˆij,u ˆi,u ˆj) and «(r ˆij,u ˆi,u ˆj) are the range and strength parameters, respectively ~see Ref. 21 for explicit expressions!, s and « depend on the anisotropy parameters k ~molecular elongation!and k 8~energetic anisotropy!. s 0and «0~the range and energy values in the side-by-side arrangement!are used to define the pressure and temperature scales in our simulations. Thus P*5P s 0 3/«0and T*5kT/«0. The dipole–dipole interaction is given by Udd~r ˆij,u ˆi,u ˆj!5@ m i m j23~ m i•r ˆij!~ m j•r ˆij!# rij 3,~3! where m i5 m u ˆiand m j5 m u ˆjdenote the dipole moments of molecules iand j, respectively. The dipoles are longitudinal, i.e., parallel to the unit vectors which represent the principal molecular axes u ˆior u ˆj. While the GB potential UGB(r ˆij,u ˆi,u ˆj) can be truncated and shifted at some distance, rc, less than L/2, Lbeing the length of the simulation box, the same is not true for the dipolar interaction. The long-ranged nature of the dipolar potential is taken into account with the reaction field method. This scheme22 assumes that particles beyond a cutoff distance, rc, act as a dielectric continuum of dielectric constant «RF .23 The interaction energy for a pair of dipoles can then be written U mm ~r ˆij,u ˆi,u ˆj!5Udd~r ˆij,u ˆi,u ˆj!22~«RF21! 2«RF11 m i• m j rc 3, ~4! for rij,rcand U mm (r ˆij,u ˆi,u ˆj)50 for rij.rc. III. SIMULATION DETAILS The simulations were performed using Monte Carlo in the isothermal–isobaric ensemble ~i.e., constant number of molecules N, pressure P, and temperature T!with the reaction field method. Each Monte Carlo cycle consists of Ntrial displacements and reorientations of the GB molecules and approximately one trial volume change. The maximum displacement, rotation, and volume change are altered to ensure that around 40% of each type of move is accepted; this value should lead to the most efficient sampling of phase space. The reaction field is a simple but accurate method for accounting, in an average way, for the long range of the dipolar interaction. The technique is much faster than the Ewald summation method, but yields results which are essentially indistinguishable. For comparison of these two methods see Refs. 9, 24, and 25. In order to analyze the orientational order and the possible polarization of the mesophases, we have calculated two order parameters Sand P1; the former ~the nematic order parameter!is obtained as the largest positive eigenvalue of the Qtensor,26 the elements of which are defined as Q ab 51 N( i51 N 1 2~3u a iu b i2 d ab !~5! ~uk iis the kth component of the axial orientation u ˆiof the ith molecule!, with the nematic director nbeing the corresponding eigenvector. The polarity P1is obtained as P151 N U ( i51 N u ˆi•n U .~6! The nematic order parameter, S, is zero in the isotropic phase, and has a value of one in a perfectly aligned phase; it provides no information as to positional order within the system. In real systems, Sattains values around 0.4 at the I–N transition, whereas in simulation, higher values are frequently encountered due to finite size effects. The polarity P1is also zero in the isotropic phase, can only attain its maximum value ~one!in a perfectly aligned state, and similarly gives no information on the positional order within the system. The functions are distinct because P1lacks ‘‘head– tail’’ symmetry, and is only nonzero when more of the dipoles ‘‘point’’ in one direction than another. As such, P1 measures the spontaneous polarization of the phase. We calculate various pair distribution functions, e.g., the first-rank orientational correlation function defined by g1~r!5 ^ P1~u ˆ1,u ˆ2! & 5 ^ cos u 12 & ,~7! and the second-rank orientational correlation function, g2~r!5 ^ P2~u ˆ1,u ˆ2! & 5 K 3 2cos2 u 1221 2 L ,~8! where u 12 is the angle between the principal molecular axes of molecules 1 and 2. In addition, we also determine the orientationally averaged pair correlation functions for two molecules whose centers lie on a line parallel, g i (r i ), or on a line perpendicular, g'(r'), to the director. Here the distances r i and r'are measured parallel and perpendicular to the director, respectively. 9531J. Chem. Phys., Vol. 109, No. 21, 1 December 1998 Houssa, Rull, and McGrother Downloaded 24 Sep 2003 to 150.214.138.210. Redistribution subject to AIP license or copyright, see http://ojps.aip.org/jcpo/jcpcr.jsp These functions alone cannot unequivocally determine the nature of certain smectic phases. To verify the structure of these liquid crystals, we calculate the bond-orientational order within the smectic layers. This function is defined as B51 6N K ( j51 6 exp~6i f ij! L ,~9! where f ij is the angle between the bond linking particles i and jand a fixed reference axis. The sum is over only nearest-neighbor bonds. We define nearest-neighbor bonds to be those within rb'1.2 s 0. This order parameter takes values close to zero when no in-plane bond order exists and close to one in the presence of perfect hexatic in-plane bond order. Principally this parameter is used to distinguish the smectic-A phase from the more ordered smectic-B structure. For GB molecules with k 53, we have simulated N 5256 molecules with longitudinal point dipoles located at the center of the molecules, with dipole moment m *. In order to minimize the effect of the simulation cell on the observed phases, a cubic simulation box was employed with periodic boundary conditions. For more elongated particles ( k 54) we are forced to examine larger systems (N5500), in a cuboidal cell. Dipole moments ranging from m *50.5 to 2.5 ~in steps of 0.5!were considered. For the smaller system, at each value of m *, simulations were initiated from a facecentered-cubic ~fcc!lattice with five layers perpendicular to the z-axis, which is expanded to a low density r *5 rs 0 3 50.08 (P*50.5) and quickly loses positional and orientational order. The system was then slowly compressed in reduced pressure steps of 0.5 or less. For the larger system, states in the isotropic, nematic, and smectic-A phase were generated with nonpolar GB systems. The dipole is introduced and the system allowed to equilibrate once more. Typically, between 3 and 53105cycles were performed for each state point, increasing to 106cycles in the vicinity of phase transitions. By modern standards these system sizes are moderate; however, we performed several simulations with larger system sizes and noted no systematic difference in the results. Furthermore, simulating small systems allows us to explore phase space much more thoroughly. At the liquid crystalline phase transitions, large molecular reorganizations must occur and lengthy runs are crucial. We are confident that system sizes are sufficient to determine the nature of the phases and FIG. 1. Phase diagram at T*51.25, of a system of 256 Gay–Berne particles ~ k 53, k 855!with embedded central, longitudinal point dipoles. The reduced pressure P*5P s 0 3/«0is plotted as a function of the reduced number density r *5 rs 0 3. Each plot depicts a different reduced dipole moment m * 5( m 2/«0 s 0 3)1/2 ~a! m *50.5, ~b! m *51.0, ~c! m *51.5, and ~d! m *52.0. Error bars denote one standard error in the density. Horizontal lines denote the position of the phase transitions. 9532 J. Chem. Phys., Vol. 109, No. 21, 1 December 1998 Houssa, Rull, and McGrother Downloaded 24 Sep 2003 to 150.214.138.210. Redistribution subject to AIP license or copyright, see http://ojps.aip.org/jcpo/jcpcr.jsp accurately predict the thermodynamic properties of the system. This belief is given credence by the excellent agreement noted between our simulations and those of Brown et al.28 In this paper we present the results of a series of NPT Monte Carlo ~MC!simulations on a system of molecules interacting via the Gay–Berne dipolar potential with k 53or 4 and k 855. IV. SIMULATIONS RESULTS A. k 53 1. T *51.25 At this temperature we have constructed isotherms for various dipolar strengths. For each value of m *, the lowest pressure simulated corresponds to an isotropic fluid, and increasing the pressure eventually leads to at least one phase transition, indicated by a small discontinuity in the density and a marked increase in orientational order. Focusing our attention on weaker dipole moments m *<2.0, the numerical phase diagrams are presented in Fig. 1, with the pressure dependence of the orientational order parameter Sand polarity P1given in Fig. 2. From these figures we can perceive two distinct transitions. The first transition is from Svalues around zero ~isotropic phase!up to nonzero values typical of the nematic phase. That the resultant phase is nematic, can be confirmed by inspection of the projections of the pair distribution functions parallel g i (r i ) and perpendicular g'(r') to the director ~Figs. 3 and 4!. These figures show that the structure is liquidlike and that there is no significant layering at the lower of the two indicated pressures. The second transition is from moderate values of S, up to saturation S'0.94. This high density phase is best characterized by g i (r i )~Fig. 3!: the sinusoidal variation of this function indicates that there is one-dimensional layering of particles, which is the trademark of the smectic liquid crystals. The wavelength of this periodic function is indicative of the layer spacing of the smectic strata and the amplitude is related to how well defined the layers are. It is clear that the strength of the dipole does not significantly affect either the wavelength or the amplitude of this function at the nematic–smectic transition. This result is quite surprising; previous simulation studies have indicated the significant impact of the dipole upon the details of the structure, particularly of the ordered phases. The striking similarity of these figures shows that the structure of the smectic liquid crystal is dependent principally on the GB interactions. The layer spacing ~see Fig. 3!is slightly less than one molecular length, which is a common feature of such phases when the constituent molecules have an ellipsoidal core: the tapering of the molecule permits a degree of interdigitation and this FIG. 2. Variation of the nematic order parameter S~triangles!and polarity P1~circles!, as a function of density for the dipolar Gay–Berne fluid. The four figures correspond to the same dipole moments as Fig. 1. Error bars denote one standard error in the order parameters. 9533J. Chem. Phys., Vol. 109, No. 21, 1 December 1998 Houssa, Rull, and McGrother Downloaded 24 Sep 2003 to 150.214.138.210. Redistribution subject to AIP license or copyright, see http://ojps.aip.org/jcpo/jcpcr.jsp structure has the most efficient packing of particles. From Fig. 4 the N–Sm transition is also accompanied by an increase in the perpendicular order. Significantly, a peak develops at less than one diameter, corresponding to nearest neighbors in the next smectic stratum. The long-ranged nature of this function in the smectic phase indicates that the phase is more highly structured than a simple smectic-A phase. In Fig. 5 we show a snapshot of the m *51.5 system at a reduced pressure P*59.0, in this smectic phase. The smectic layers are clear from Fig. 5~a!, and the positional order within those layers is evident from Fig. 5~b!. The N–Sm transition is accompanied by a distinct increase in the bond-orientational order. Thus we conclude that the phase is smectic-B. Interestingly, as the strength of the dipole is increased, the number of peaks perpendicular to the director is seen to diminish. The orientationally averaged pair distribution function, Fig. 6, exhibits the expected oscillatory behavior and high value at long distances, of the layered smectic phases. From Fig. 6 we see that there is significant shortrange orientational order in the isotropic phase near to the I–N transition. This order is seen to increase slightly with increasing dipole moment. After the I–N transition these orientational correlations no longer fade to zero within the simulation cell. The orientational correlation at contact is an increasing function of m . In Fig. 7 we display the first-rank orientational correlation function for states in the isotropic nematic and smectic-B phases for these four dipole moments. Negative values indicate that particles are antiparallel. The figures are qualitatively equivalent; each shows that: nearest neighbors are antiparallel; next nearest neighbors are oriented parallel to the selected particle; the minima and maxima increase in magnitude as one moves from isotropic to nematic to smectic phases; the order dies out within the simulation cell, indicating no globally polar phases. Both of the transitions ~I–N and N–SmB!are seen to be weakly first order. Neither phase change affects the polarity P1, which remains essentially zero, over the entire pressure range. Overall, these systems have the phase sequence I–N– SmB on compressing from low density at this temperature (T*51.25). This is true for nonpolar GB particles and those with weak, central, longitudinal dipoles. A natural question is why the saturation value of Sis less than one. A perfectly aligned state is never achieved in simulations as a result of several particles becoming trapped at right angles to the director. These particles are particularly noticeable in the smectic phases, where they prefer to position themselves in the interlayer region. This fact is the source of the characteristic minimum in the orientational pair distribution function of simulated smectic liquid crystals. FIG. 3. The projection of the radial pair distribution function parallel to the director g i (r i ) for states on either side of the nematic–smectic transition for the same systems as Fig. 1. 9534 J. Chem. Phys., Vol. 109, No. 21, 1 December 1998 Houssa, Rull, and McGrother Downloaded 24 Sep 2003 to 150.214.138.210. Redistribution subject to AIP license or copyright, see http://ojps.aip.org/jcpo/jcpcr.jsp The behavior of these molecules has been studied in detail for related systems.27 In Table I we report our estimates of the I–N and N– SmB transition pressures and coexisting densities for the different dipolar strengths. From this table, it is clear that the introduction of the longitudinal dipole into the model has some impact on the I–N transition pressure, but the strength of the dipole has little effect upon the density of the isotropic phase at the transition. The density of the nematic phase that coexists with the isotropic fluid is also unchanging with variations in dipole moment. The I–N transition pressure (PIN) is seen to fall with increasing dipole moment. The consistency of the I–N transition densities r IN , which are, to within experimental error, identical for each of the dipole strengths which lead to an I–N transition ( m *<2.0), is striking. At the layering transition ~N–SmB!, we now see that the densities at coexistence decrease as the strength of the multipolar interaction is increased. The transition pressure shows a very clear trend, moving to lower values as the dipole moment is increased. This is in accordance with expectation: longitudinal dipoles are well known to promote layered structures in both numerical and experimental studies, thus the stabilization of the smectic phase is to be expected. Increasing the dipole makes the smectic-B phase stable at lower pressures. Configurations at P*57.0 were used as starting points for decreasing the pressure. A limited number of expansion runs have been performed to gauge the level of hysteresis at the transitions. For these weak dipoles ( m *<2.0), we see no evidence of hysteresis at the I–N transition; however, at the N–SmB phase change, the transition pressure is very different depending on the direction from which the transition is approached. The same is true of the N–SmB transition density, but since the phase diagram is very steep in this region ~i.e., large changes in pressure correspond to only small changes in density!, the effect is less obvious. The natural consequence of the trends noted for the weaker dipoles in Table I is that the nematic phase will disappear at a triple point as the strength of the dipole moment is increased. This is due to the increasing stabilization of the smectic-B phase by the dipole, leading eventually to the layered phase preempting the nematic phase when the multipole is sufficiently strong. Figure 8 is the phase diagram for m *52.5. A strongly first-order phase transition, directly from the isotropic fluid to the smectic-B phase, is observed. The nature of the transition is made clear by a jump in the bond orientational order at the I–LC transition. Alignment of the molecules is accompanied not only by layering, but also the development of long-ranged bond order within these layers. The bond-order parameter takes a value of B50.55 at the transition. For these highly polar GB molecules, the I–SmB transition oc- FIG. 4. The projection of the radial pair distribution function perpendicular to the director g'(r') for the same states as Fig. 3. 9535J. Chem. Phys., Vol. 109, No. 21, 1 December 1998 Houssa, Rull, and McGrother Downloaded 24 Sep 2003 to 150.214.138.210. Redistribution subject to AIP license or copyright, see http://ojps.aip.org/jcpo/jcpcr.jsp FIG. 5. Snapshot of the dipolar Gay–Berne fluid at T*51.25, P*59.0, m *51.5, in the smectic-B phase. The same configuration is shown from different perspectives: ~a!from the side to show the smectic layering and ~b!exhibiting the strong positional correlation between successive layers. The colors indicate the direction of the dipole. The size of the particles is reduced for clarity. FIG. 6. Second-rank orientational correlation function g2(r) for dipolar Gay–Berne fluids at the I–N and N–SmB transitions. The four figures represent different dipole moments ~as for Fig. 1!. 9536 J. Chem. Phys., Vol. 109, No. 21, 1 December 1998 Houssa, Rull, and McGrother Downloaded 24 Sep 2003 to 150.214.138.210. Redistribution subject to AIP license or copyright, see http://ojps.aip.org/jcpo/jcpcr.jsp curs at a significantly lower pressure than the I–N transition for the weaker dipolar systems. Once again, though, we note from Table I that the density at which the isotropic phase ceases to be stable is unaffected r *'0.305. The position of the phase transition is inferred from the pressure versus density phase diagram, Fig. 8. Even more compelling evidence for the location of the phase transition is provided by the behavior of the nematic order parameter S and polarity P1as a function of density ~Fig. 9!. The variation of the order parameters also provides information as to the nature of the resultant phases. The first nonzero value of Soccurs at r *'0.348, but now we can see notable pretran- FIG. 7. First-rank orientational correlation function g1(r) for dipolar Gay–Berne fluids in the isotropic, nematic, and smectic-B phases. The four figures represent different dipole moments ~as for Fig. 1!. FIG. 8. Phase diagram at T*51.25, of a system of 256 Gay–Berne particles ~ k 53, k 855!with embedded central, longitudinal point dipoles. The reduced pressure P*5P s 0 3/«0is plotted as a function of the reduced number density r *5 rs 0 3. The reduced dipole moment is m *5( m 2/«0 s 0 3)1/252.5. Error bars represent one standard error in the density. TABLE I. Transition pressures and densities for dipolar GB fluid at a reduced temperature of T*51.25 as obtained by MCNPT simulation. GB parameters are: k 53, k 855. m *is the reduced dipole moment, PIN is the isotropic to nematic transition pressure, PNSmB is the same property for the nematic to smectic-B transition, r IN are the densities of the coexisting isotropic and nematic phases, and r NSmB are the same values for the layering transition. Note that the strongest dipole has a direct isotropic to smectic-B transition. The system size is N5256. m *PIN r IN PNSmB r NSmB 0.5 4.99 0.304, 0.319 9.49 0.368, 0.380 1.0 5.00 0.305, 0.320 9.00 0.364, 0.379 1.5 5.00 0.307, 0.322 7.00 0.346, 0.361 2.0 4.49 0.305, 0.321 4.99 0.321, 0.339 PISmB r ISmB 2.5 3.70 0.300, 0.348 9537J. Chem. Phys., Vol. 109, No. 21, 1 December 1998 Houssa, Rull, and McGrother Downloaded 24 Sep 2003 to 150.214.138.210. Redistribution subject to AIP license or copyright, see http://ojps.aip.org/jcpo/jcpcr.jsp