Molecular Dynamics Simulation Study of the Influence of Chirality on

Jun 27, 2008 - Liquid Crystal Nano-system Project, ERATO/SORST, Japan Science and Technology Agency, 5-9-9 ... B , 2008, 112 (29), pp 8452–8458...
0 downloads 0 Views 2MB Size
8452

J. Phys. Chem. B 2008, 112, 8452–8458

Molecular Dynamics Simulation Study of the Influence of Chirality on the Stability of the Smectic Q Liquid Crystal Phase Makoto Yoneya,*,†,‡ Takahiro Yamamoto,†,‡ Isa Nishiyama,† and Hiroshi Yokoyama†,‡ Liquid Crystal Nano-system Project, ERATO/SORST, Japan Science and Technology Agency, 5-9-9 Tokodai, Tsukuba 300-2635, Japan, and Nanotechnology Research Institute, National Institute of AdVanced Industrial Science and Technology, 1-1-4 Umezono, Tsukuba 305-8568, Japan ReceiVed: December 28, 2007; ReVised Manuscript ReceiVed: April 11, 2008

The structure of smectic Q (SmQ) liquid crystal phase consisting of a dichiral molecule, called M7BBM7, was studied by submicrosecond molecular dynamics (MD) simulation. A detailed atomic model was used to study the stability of a model SmQ structure proposed by Levelut et al. (Levelut, A.-M.; Hallouin, E.; Bennemenn, D.; Heppke, G.; Lotzsch, D. J. Phys. II 1997, 7, 981) and its difference between (S,S)-, (S,R)M7BBM7 and racemic mixture systems. Negative values of the fourth-rank orientational order parameter (〈P4〉), which characterize the model SmQ structure, were stably kept up to a 100 ns MD run only in the (S,S)-M7BBM7 system and lost in the other systems. The results correspond well to the marked chiral sensitivity in real systems where only the (S,S)-M7BBM7 system (among the three above-mentioned systems) shows the SmQ phase. Our simulation results imply that the asymmetric intramolecular potentials and resultant chirality-dependent molecular conformations are primarily responsible for keeping the negative values of 〈P4〉 and the model SmQ structure. I. Introduction Molecular chirality is a subject of great interest in the chemical physics community and plays a significant role in various matters including liquid crystals (LCs). Blue phases (BPs), which are LC phases with three-dimensional (3D) nanostructures, are considered to appear with strong molecular chirality, which leads to local twist (helix) formation in two orthogonal directions and a further double twist cylinder as a building block of BP’s 3D structure.1 The BP unit cell has cubic symmetry with a cell length close to the helical pitch and is on the order of 102-103 nm. Recently, the BPs have attracted increasing interest with respect to their emerging applications, such as BP lasers2 and BP displays.3 Smectic Q (SmQ) phases4,5 are also considered as chiralityoriginated 3D-structured LC phases since the SmQ phases are very sensitive to optical purities and are lost with low chiral fractions.6 Only a limited number of compounds have been known to exhibit the SmQ phases, but new compounds have been reported recently.7–11 Compared with the BPs, the SmQ phases have a small unit cell length (∼10 nm) and tetragonal (for T*1, T*2, and T*3 SmQ subphases) or hexagonal (for the H* subphase) symmetry.5 The unit cell length of a SmQ phase is similar to that of yet another 3D-structured LC phase, the cubic (Cub) phase, although the origin of Cub phase formation is considered to be different, that is, nanosegregation.12 However, there is experimental evidence that implies some connections between SmQ and Cub phases. First, some SmQ and Cub compounds have been reported to be miscible as a binary mixture, and this mixture shows a SmQ-Cub phase transition.4 Second, another SmQ compound shows Cub and the tetragonal phase (called as S4 phase)13 in its racemic mixture.5 In the latter case, the space groups of the SmQ, Cub, and S4 phases have * To whom correspondence should be addressed. † Japan Science and Technology Agency. ‡ National Institute of Advanced Industrial Science and Technology.

Figure 1. Relationships between Cub(gyroid), S4, SmQ(T*2), and BP1 phases.

been reported to be I4122, Ia3d, and I41/acd, respectively. These space groups are closely related to each other, as depicted in Figure 1. I41/acd (S4) is considered to be a tetragonal modification of cubic Ia3d (double gyroid).14 I4122 (SmQ:T*2) corresponds to the chiral analogue of I41/acd and also the tetragonal modification of chiral cubic I4132 (single gyroid), which is a space group of BP1. Relationships similar to those in Figure 1 could be found between SmQ(T*1)-Cub(primitive) and SmQ(T*3)-Cub(diamond) (see the Supporting Information).15 Clarifying the relationships between the SmQ and Cub phases would be a key for understanding the nanoscale 3D structure formation in LC systems. We have been trying the molecular dynamics (MD) simulation approach to investigate the molecular level structure of the Cub phases.16,17 If we could reasonably simulate both the SmQ and Cub phase structures, we could further analyze their precise relationships. From the viewpoint of molecular simulation of liquid crystals, chirality-related phenomena are challenging subjects even in the simplest cholesteric (chiral nematic) phase.18–20 One reason for this is that the chiral interactions are considered to be subtle, and the resultant structures, that is, cholesteric helices, are generally on much larger scales (∼few microns) than the scale that the molecular simulation can currently handle (efew tens

10.1021/jp712121v CCC: $40.75  2008 American Chemical Society Published on Web 06/27/2008

Influence of Chirality on the Stability of the SmQ

J. Phys. Chem. B, Vol. 112, No. 29, 2008 8453

Figure 2. (S,S)-M7BBM7 (Cry - 329 K - SmC*A - 333 K - SmQ(T*1) - 359 K - Iso).

Figure 3. Fitted torsional potentials (line) to the ab initio MO results ( ×) for rotation around the bond displayed in the inset.

of nanometers).20 However, in the case of the SmQ phase, its strong chiral sensitivity and ∼10 nm order unit cell scale may make the study of the chiral effects by molecular simulation possible. In this study, the SmQ phase structure of a dichiral compound biphenyl-4,4′-dicarboxylic acid 4′-(1-methylheptyl)ester 4-[4′(1-methylheptyloxycarbonyl)-biphenyl-4-yl]ester (called M7BBM7;6 shown in Figures 1 and 2) was studied by MD simulation. As a first step, we report the result of the MD simulation study on the stability of a previously proposed SmQ structural model5 and its dependence on the chirality of a system, for example, a pure enantiomer or racemic mixture. II. Model and Method The molecular model utilized was a detailed atomic model using “united” atoms, that is, all-hydrogen atoms in CHn groups were treated as “united” atoms, except the hydrogen atoms attached to the chiral carbon atoms.21 Each molecule was treated as a flexible model, except that all of the bond-stretching degrees of freedom were constrained to its equilibrium bond length. The inter- and intramolecular interaction potential functions employed were basically GROMACS potentials.22 The torsional potentials around the -C*-O- bond in the chiral part of the M7BBM7 molecule were fitted to ab initio molecular orbital (MO) calculation (with HF/6-31G* level) results utilizing the program Gaussian98, as shown in Figure 3. The primary fitting (including chiral asymmetry) was performed by adjusting 1-4 Lennard-Jones intramolecular interaction parameters related to this torsional degree of freedom, and then, the residual (symmetric) component was fitted with the following torsional angle function V(φ) 5

V(φ) )

∑ Cn cosn(φ)

n)0

The above-mentioned modeling enables a change in molecular chirality simply by changing the geometrical coordinates of the

Figure 4. DTCBCL structure representations with eight unit clusters (a) and four anticlinic bilayer clusters (b). The molecules drawn with “nail heads” are tilted out of the page, with the nail head ends toward the front of the page.

atoms attached to the chiral carbon atom. A similar fitting was performed for the inter-ring torsional potential of biphenyl parts. Atomic charges were also obtained from the Gaussian98 calculation results. The GROMACS molecular topology file for the M7BBM7 molecule with atomic charge assignments that help reproduce the simulation results in this study is found in the Supporting Information. Trajectories were produced by the MD program package GROMACS (version 3)22 with the leapfrog time integration and LINCS bond constraint,23 and the time integration step was set to be 2 fs. Charge-group-based twin-range 0.9 nm van der Waals and 1.8 nm electrostatic cutoff distances22 were applied to nonbonded interactions. The initial structure for MD simulation was constructed from the model SmQ(T*1) structure proposed by Levelut et al.,5 which is called “double-twist clusters on a body-centered lattice (DTCBCL)” (Figure 4). The original DTCBCL model (Figure 4a) can be represented with four anticlinic bilayer clusters (Figure 4b). We note that this SmQ model is similar to our gyroid Cub model16 in its structural hierarchy with the clusters as the basic units. To obtain this anticlinic bilayer cluster, first, we performed the MD simulation at the SmC*A phase temperature (330 K) with 256 molecules (1:1 mixture of the (S,S)- and (R,R)-M7BBM7 models). The snapshot image of the simulated SmC*A-like structure after a 32 ns MD run is shown in Figure 5. From this SmC*A-like structure, 126 M7BBM7 molecules at its MD cell center region were extracted to form the unit anticlinic bilayer cluster in Figure 4b, and then, the entire DTCBCL structure (with 126 × 4 ) 504 molecules) was constructed, as shown in Figure 6.

8454 J. Phys. Chem. B, Vol. 112, No. 29, 2008

Yoneya et al.

Figure 5. Simulated SmC*A-like structure.

The number of molecules (504) corresponds to the value deduced from the corresponding SmQ phase cell parameter (a ) b ) 8.04 nm, c ) 8.64 nm) and the estimated density of 1.0 g/cm3. To examine the effect of chirality on the stability of the model SmQ(T*1) structure, three different systems were constructed from the system in Figure 6 by changing the chiralities of the M7BBM7 model molecules. First, a pure (S,S)-M7BBM7 system was obtained by changing all of the rest of the half of the (R,R)-M7BBM7 molecules to the (S,S)-M7BBM7 molecules. Second, a pure (S,R)-M7BBM7 system was obtained by changing either side of two chiral parts of all of the M7BBM7 model molecules. Third, a racemic mixture system with equal numbers of the four isomers,24 that is, (S,S)-, (R,R)-, (S,R)-, and (R,S)-M7BBM7 molecules, was obtained by changing the chirality of the (randomly selected) M7BBM7 model molecules. Among these three systems, only the (S,S)-M7BBM7 system shows a SmQ phase in the real compound, and the other (S,R)M7BBM7 and racemic mixture systems do not.6,25 The reason why the (S,R)-M7BBM7 system does not show a SmQ phase was considered to be its very low chirality (2 orders of magnitude smaller helical twisting power and spontaneous polarization values than those of the (S,S)-M7BBM7 system).25 SmQ temperature simulations were performed, starting from these systems placed in a tetragonal simulation box under periodic boundary conditions. After a few thousand steps of the steepest descent energy minimization and (∼4 ns) preproduction MD runs, production MD runs up to 100 ns were carried out at a constant average temperature and a pressure under weak couplings to the heat and pressure bath26 (with time constants of 0.1 and 4.0 ps, respectively). The simulation temperature was set at around the middle of the SmQ phase temperature range of the (S,S)-M7BBM7 system, that is, 350 K (see the caption of Figure 2). This corresponds to the SmC*A temperature range for the (S,R)-M7BBM7 and racemic mixture systems.6,25 Simulation pressures were set to ensure the average simulated densities at around 1.0 g/cm3 for the (S,S)-M7BBM7 system, and the same settings were applied to all three above-mentioned systems. The above-mentioned pressure setting was actually negative, that is, -28.0 MPa, as that in our previous work.16 This negative pressure setting was required due to the known deficiency in GROMACS (GROMOS) force field.27 The computational time was about 17 h per 1 ns (5 × 105 time steps) of simulation on the 64 node Pentium III/1 GHz cluster connected via Myrinet2000. III. Results Temporal variations in densities in the MD runs at the SmQ temperature (350 K) starting from the DTCBCL structure

Figure 6. Initial DTCBCL model structure with four anticlinic bilayers [(a) simple stacking and (b) after applying periodic boundary conditions]. Only mesogenic core parts of each molecules are shown for clarity.

(Figure 6) with the (S,S)-, (S,R)-M7BBM7, and racemic mixture systems are shown in Figure 7. The simulated densities appear to be well equilibrated at around the end of the MD runs (100 ns) for all of the systems, but they differ depending on the molecular chiralities of the systems. The (S,R)-M7BBM7 system keeps a slightly larger density than the other systems. In contrast, the racemic mixture system initially appears to be equilibrated at around slightly lower densities than the (S,S)-M7BBM7 system, but they become equilibrated at similar densities after a 60 ns MD run. The latter change in the racemic mixture system may imply a change in molecular arrangement to homochiral-like local packing. A. Orientational Structures. To analyze the orientational structure of a simulated system, we evaluated the second- and fourth-rank orientational order parameters denoted as 〈P2〉 and 〈P4〉, respectively. If we assume the ideal DTCBCL structure

Influence of Chirality on the Stability of the SmQ

J. Phys. Chem. B, Vol. 112, No. 29, 2008 8455

Figure 7. Time evolutions of the densities of the pure (S,S) and (S,R) enantiomers and racemic mixture systems in MD runs from the DTCBCL initial structure.

Figure 8. Molecular tilt angle θ dependencies of the second- and fourth-rank orientational order parameters of the DTCBCL structure.

(Figure 4), the corresponding 〈P2〉 and 〈P4〉 can be calculated as functions of the molecular tilt θ simply from the Legendre polynomials P2(θ) ) (3 cos2θ - 1)/2 and P4(θ) ) (35 cos4θ 30 cos2θ + 3)/8, respectively. In this case, the molecular tilt angle θ is defined as the angle with respect to the c-axis in Figure 4 (normal to the page). The obtained P2(θ) and P4(θ) are plotted in Figure 8. The tilt angle of around 45° was experimentally reported for the (S,S)-M7BBM7 compound.28 With this tilt angle in Figure 8, 〈P2〉 and 〈P4〉 would be small (positive) and negative, respectively. Then, with reference to 〈P4〉, we can distinguish the DTCBCL structure (gives 〈P4〉 < 0) from the spherical (isotropic) orientational symmetry structure (gives 〈P4〉 ∼ 0). 〈P2〉 was calculated from the largest positive eigenvalue of the order parameter tensor30 Q

QRβ )

1 2N

∑ 3eiReiβ - δRβ i

which consists of molecular axis vectors i that are determined by the connecting vectors between two carbon atoms at each end of the outermost phenyl rings. By using the eigenvector that corresponds to the above-mentioned eigenvalue, 〈P4〉 was calculated from the fourth-rank tensor, as in our previous work.16

Figure 9. Time evolutions of (a) the second- and (b) fourth-rank orientational order parameters in MD runs from the DTCBCL initial structure.

However, the common selection of the largest positive eigenvalue as 〈P2〉30 fails in the case of θ > 54.7° (called the magic angle) in Figure 8 as P2(θ) becomes negative. In contrast, the selection by Romano,29 that is, the largest absolute eigenvalue, can resolve this problem and is valid even for θ > 54.7°. We had tried both selections to check their difference and found no significant difference between them for all three simulated systems. This implies that the simulated tilt angles did not largely exceed the magic angle. The calculated 〈P2〉 and 〈P4〉 values are plotted in Figure 9. The second-rank order parameter 〈P2〉 showed a very small value for all of the (S,S)-, (S,R)-M7BBM7, and racemic mixture systems because of the high tilt angles (∼50°) within the simulated unit anticlinic bilayer. In the case of the fourthrank order parameter 〈P4〉, the values were negatively shifted in the (S,S)-M7BBM7 system, whereas these were around zero in the other systems. This indicates that the negative values of the fourth-rank orientational order parameter that characterize the DTCBCL structure were stably kept only in the (S,S)-M7BBM7 system. The notable stability of the DTCBCL model characteristics with the (S,S)-M7BBM7 system is not trivial because the slightly different model

8456 J. Phys. Chem. B, Vol. 112, No. 29, 2008

Yoneya et al.

Figure 10. Structure in which the molecular tilting directions with respect to the c-axis are different from those of the DTCBCL model, that is, the molecules are lying flat on the page.

Figure 12. Mean-squared displacements in the last 12 ns (88-100 ns) MD run for the (S,S)-M7BBM7 system (upper: total; lower: for x, y, and z directions).

Figure 11. Comparison of the structure functions (enlarged in the inset) for the initial DTCBCL structure and simulated (S,S)- and (S,R)M7BBM7 structures averaged over the last 100 ps period of 100 ns MD runs.

structure (in Figure 10) was not stable with the same simulation model and procedures. We constructed the structure (with the (S,S)-M7BBM7 molecules) depicted in Figure 10 in which the molecular tilting directions with respect to the c-axis are different from those in Figure 4 (these two different tilting manners have both been discussed for the SmQ phase structures).5 Corresponding to the MD run in the previous section starting from the structure in Figure 4, an equivalent run starting from the structure in Figure 10 was performed for comparison. We found that the latter structure transformed to a uniformly oriented SmA-like structure within the 1 ns MD run, that is, only the MD run starting from the former structure (Figure 4), was stable up to 100 ns. This result suggests that the structure in Figure 4 is the more likely the real structure of the SmQ phase. Actually, the experimental study, where the tilting manner corresponding to that in Figure 4 would be more probable than that in Figure 10, was reported in a particular SmQ system.8 B. Positional Structures. Next, we calculated the structure function S(q) from the simulated structures to analyze the positional structures, where q is the scattering wave vector. Figure 11 shows a detailed comparison of the averaged structure functions between the initial model structure and the last simulated structure of the three systems.

As in this figure, the peaks of the initial structure function are well (and best among the three) preserved in the (S,S)M7BBM7 system even after 100 ns. Peak shifts compared to the initial structure, most evident at around q ∼ 5 nm-1 from the initial structure may not be directly related to the changes in the density in Figure 7 since all three systems show similar shifts. One possible cause would be the structural relaxation from the constructed initial structure. Molecular mean square displacements from simulated trajectories in the last 12 ns are plotted in Figure 12 for the (S,S)M7BBM7 system. The obtained plot shows a steady molecular diffusion, indicating that the simulated state is not a crystal or glassy state but a liquid crystalline state, although the state has a 3D structure, as in the above-mentioned structure functions. The difference between the separated plots for the x, y, and z directions is very small, that is, the diffusion is almost isotropic. The diffusion constant obtained by fitting the line (in the last 10 ns period) is 0.06 × 10-11 m2/s. This is 1 order of magnitude smaller than that from the NMR study31 (3.0 × 10-11 m2/s) for the thermotropic cubic mesogen, 4′-n-hexadecyloxy-3′-nitrobiphenylcarboxylic acid (ANBC16). IV. Discussion As in the Introduction, SmQ phases are very sensitive to chirality and are lost in the racemic mixture and low-chirality compounds (e.g., (S,R)-M7BBM7). Our simulation results show a marked difference between the (S,S)-M7BBM7 system and the racemic mixture and (S,R)-M7BBM7 systems. One possible origin of this marked difference could be chirality-dependent molecular conformations. The large asymmetry of the torsional potentials at the M7BBM7 chiral part (Figure 3) may cause specific chirality-dependent conformations. Actually, the distributions of the torsional angles near this chiral part differ depending on the molecular chiralities of the systems as in Figure 13; thus, the corresponding molecular conformations also differ depending on the chiralities. The density (packing) differences in Figure 7 could be related to the chirality-dependent molecular conformations. The characteristics of the chiral part of the M7BBM7 molecule would be important since all known SmQ-phase-showing compounds commonly have the same chiral part structure.4–11 In our study, the negative values of the fourth-rank orientational order parameter 〈P4〉 that characterize the model

Influence of Chirality on the Stability of the SmQ

J. Phys. Chem. B, Vol. 112, No. 29, 2008 8457 the chiral part. The time evolution of 〈P4〉 in Figure 14 reveals that an initially positive 〈P4〉 became negative with the intra(S,S)-M7BBM7 model, whereas it remained positive with the inter-(S,S)-M7BBM7 model within a 3 ns MD run. The results imply that the asymmetric intramolecular potentials and resultant chirality-dependent molecular conformations are primarily responsible for attaining a negative 〈P4〉 and the model SmQ (DTCBCL) structure. In this paper, we describe our molecular simulation results regarding the stability difference of the model SmQ structure with three different chirality systems. The results show a positive prospect in the study of complex chiral effects with the atomistic molecular simulations, although we are only at the starting point.

Figure 13. Torsional angle distributions (averaged over 100 ps in t ) 1.9-2.0 ns) for three systems with different molecular chiralities.

Acknowledgment. We thank the BIRD Division of the Japanese Science and Technology Agency for providing computer time on their cluster computing system. This work was supported by KAKENHI (Grant-in-Aid for Scientific Research) on Priority Area “Soft Matter Physics” from the Ministry of Education, Culture, Sports, Science and Technology of Japan. Supporting Information Available: Molecular topology and additional space group relationships. This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes

Figure 14. Time evolutions of the fourth-rank orientational order parameters in the preproduction MD runs. The origin of the time (t ) 0) was adjusted to that in Figure 9b.

SmQ (DTCBCL) structure were stably kept only in the (S,S)M7BBM7 system and lost in the other systems (e.g., (S,R)M7BBM7 system). The difference in 〈P4〉 was already observed at the beginning of our MD runs (t ) 0 in Figure 9b). Actually, this difference was established in the very early stage of the preproduction (before t ) 0) MD runs, as shown in Figure 14. This rapid discrimination implies that a short time scale local change, such as a molecular conformational change, may be responsible for the discrimination. To confirm this, we performed additional MD runs starting from the configuration of the (S,R)-M7BBM7 system at t ) -3.0 ns in Figure 14 and by applying two different modified (S,S)M7BBM7 molecular models denoted as intra- and inter-(S,S)M7BBM7 models in Figure 14. In the intra-(S,S)-M7BBM7 model, the asymmetry of the intramolecular torsional potentials at the M7BBM7 chiral part (Figure 3) was basically kept, but the intermolecular asymmetry around the chiral part was lost by replacing the CH3 atoms attached to the chiral carbon atoms with hydrogen atoms. In contrast, the asymmetry of the torsional potentials was lost by replacing the torsional potentials and 1-4 intramolecular interactions with those of achiral alkyl chains in the inter-(S,S)-M7BBM7 model while keeping the intermolecular asymmetry around

(1) Dubois-Violette, E.; Pansu, B. Mol. Cryst. Liq. Cryst. 1988, 165, 151. (2) Cao, W.; Munoz, A.; Palffy-Muhoray, P.; Taheri, B. Nat. Mater. 2002, 2, 111. (3) Hisakado, Y.; Kikuchi, H.; Nagamura, T.; Kajiyama, T. AdV. Mater. 2005, 17, 96. (4) Levelut, A. M.; Germain, C.; Keller, P.; Liebert, L.; Billard, J. J. Phys. (Paris) 1983, 44, 623. (5) Levelut, A.-M.; Hallouin, E.; Bennemann, D.; Heppke, G.; Lotzsch, D. J. Phys. II 1997, 7, 981. (6) Bennemann, D.; Heppke, G.; Levelut, A. M.; Lo¨tzsch, D. Mol. Cryst. Liq. Cryst. 1995, 260, 351. (7) Nguyen, H. T.; Ismaili, M.; Isaert, N.; Achard, M. F. J. Mater. Chem. 2004, 14, 1560. (8) Pansu, B.; Nastishin, Y.; Impe´ror-Clerc, M.; Veber, M.; Nguyen, H. T. Eur. Phys. J. E 2004, 15, 225. (9) Nishiyama, I.; Yamamoto, J.; Goodby, J. W.; Yokoyama, H. Chem. Mater. 2004, 16, 3212. (10) Manai, M.; Gharbi, A.; Marcerou, J. P.; Nguyen, H. T.; Rouillon, J. C. Physica B 2005, 368, 168. (11) Yamamoto, T.; Nishiyama, I.; Yokoyama, H. Chem. Lett. 2007, 36, 1108. (12) Yoneya, M.; Araya, K.; Nishikawa, E.; Yokoyama, H. J. Phys. Chem. B 2004, 108, 8099. (13) Levelut, A.-M.; Clerc, M. Liq. Cryst. 1998, 24, 105. (14) Fogden, A.; Hyde, S. T. Eur. Phys. J. B 1999, 7, 91. (15) Johnson, C. K.; Bernett, M. N.; Dunbar, W. D. In Crystallographic Computing 7; Bourne, P. E., Watenpaugh, K., Eds.; Oxford University Press: London, 2004. (16) Yoneya, M.; Nishikawa, E.; Yokoyama, H. J. Chem. Phys. 2004, 120, 3699. (17) Yoneya, M.; Nishikawa, E.; Yokoyama, H. J. Chem. Phys. 2004, 121, 7520. (18) Allen, M. P.; Masters, A. J. Mol. Phys. 1993, 79, 277. (19) Memmer, R.; Kuball, H. G.; Scho¨nhofer, A. Liq. Cryst. 1993, 15, 345. (20) Yoneya, M.; Berendsen, H. J. C. J. Phys. Soc. Jpn. 1994, 63, 1025. (21) Yoneya, M.; Nishiyama, I.; Yokoyama, H. J. Phys. Soc. Jpn. 2003, 72, 1403. (22) Lindahl, E.; Hess, B.; van der Spoel, D. J. Mol. Model. 2001, 7, 306. (23) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. J. Comput. Chem. 1997, 18, 1463. (24) Barretto, G. M.; Collings, P. J.; Bennemann, D.; Lotzsch, D.; Hepke, G. Liq. Cryst. 2001, 28, 629.

8458 J. Phys. Chem. B, Vol. 112, No. 29, 2008 (25) Bennemann, D.; Heppke, G.; Lo¨tzsch, D.; Paus, S. In Proceedings of the 24th Freiburger Arbeitstagung; Freiburger Arbeitstagung: Freiburg, Germany, 1995; p 27. (26) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684. (27) Schuler, L. D.; Daura, X.; van Gunsteren, W. F. J. Comput. Chem. 2001, 22, 1205. (28) Lagerwall, J.; Giesselmann, F.; Seibmann, C.; Rauch, S.; Heppke, G. J. Chem. Phys. 2005, 122, 144906.

Yoneya et al. (29) Romano, S. In Physics of Liquid Crystalline Materials; Khoo, I.-C., Simoni F., Eds.; Gordon and Breach: Philadelphia, PA, 1991; Chapter 3. (30) Allen M. P.; Tildesley, D. J. In Computer Simulation of Liquids; Oxford Science Publications: Oxford, U.K., 1987; Chapter 3. (31) Ukleja, P.; Siatkowski, R. E.; Neubert, M. E. Phys. ReV. A 1988, 38, 4815.

JP712121V