- Research article
- Open Access

# Core charge distribution and self assembly of columnar phases: the case of triphenylenes and azatriphenylenes

- Silvia Orlandi
^{1}, - Luca Muccioli
^{1}, - Matteo Ricci
^{1}, - Roberto Berardi
^{1}and - Claudio Zannoni
^{1}Email author

**Received:**15 March 2007**Accepted:**08 June 2007**Published:**08 June 2007

## Abstract

### Background

The relation betweeen the structure of discotic molecules and columnar properties, a crucial point for the realization of new advanced materials, is still largely unknown. A paradigmatic case is that hexa-alkyl-thio substituted triphenylenes present mesogenic behavior while the corresponding azatriphenylenes, similar in shape and chemical structure, but with a different core charge distribution, do not form any liquid crystalline mesophase. This study is aimed at investigating, with the help of computer simulations techniques, the effects on phase behaviour of changes of the charge distribution in the discotic core.

### Results

We described the shape and the pair, dispersive and electrostatic, interactions of hexa alkyl triphenylenes by uniaxial Gay-Berne discs with embedded point charges. Gay-Berne parameters were deduced by fitting the dispersive energies obtained from an atomistic molecular dynamics simulation of a small sample of hexa-octyl-thio triphenylene molecules in columnar phase, while a genetic algorithm was used to get a minimal set of point charges that properly reproduces the ab anitio electrostatic potential. We performed Monte Carlo simulations of three molecular models: the pure Gay-Berne disc, used as a reference, the Gay-Berne disc with hexa-thio triphenylene point charges, the Gay-Berne disc with hexa-thio azatriphenylene point charges. The phase diagram of the pure model evidences a rich polymorphism, with isotropic, columnar and crystalline phases at low pressure, and the appearance of nematic phase at higher pressure.

### Conclusion

We found that the intermolecular electrostatic potential among the cores is fundamental in sta-bilizing/destabilizing columnar phases; in particular the triphenylene charge distribution stabilizes the columnar structure, while the azatriphenylene distribution suppresses its formation in favor of the nematic phase. We believe the present model could be successfully employed as the basis for coarse-grained level simulations of a wider class of triphenylene derivatives.

## Keywords

- Monte Carlo
- Radial Distribution Function
- Triphenylene
- Nematic Phase
- Molecular Electrostatic Potential

## Background

Given their capability to provide microscopic details that cannot or are difficult to access with experimental techniques, computer simulations have already demonstrated their usefulness in the field of liquid crystals [8, 9], even if the majority of studies has been carried out on calamitic-type mesogens. As liquid crystalline properties depend on collective molecular organization, the model has to contain on one hand enough atomic details to be chemically meaningful [10, 11] while on the other, still be manageable in the simulation of large enough samples of different phases [12]. When the latter aspect is predominant, the choice of simplified models which describe the molecule with a single or a few interacting sites is to be preferred. In this specific niche, simulations based on the attractive-repulsive Gay-Berne (GB) potential [13] have been very successful in reproducing and even predicting a variety of liquid crystalline phases [14, 15].

So far, however, GB simulations of discotics have not considered the real dimensions of the typical triphenylenic molecules giving columnar phases, and the aspect ratio of the actual mesogenic molecules (e.g. the core thickness/diameter ratio) has been typically overestimated. For instance, most studies [16–21], based on the parameterization of Emerson *et al*. [16] have used an aspect ratio of 0.345, taking into account only the dimensions of the rigid molecular core, disregarding the presence of alkyl chains, which actually play a crucial role in determining the mesogenic properties of these compounds. X-ray experimental studies have demonstrated [3] that a more meaningful aspect ratio, estimated as the ratio between the cofacial and the intercolumnar distances in non-tilted phases, is hardly higher than 0.2 for triphenylenes, and in some cases, like hexabenzocoronene derivatives, even lower [22]. Despite this, until now only one simulation work, by Ryckaert and coworkers, reports the phase diagram of GB discs with lower aspect ratio [20].

Considering the necessity of reparameterizing the GB potential for triphenylenes, we attempted the investigation of the differences between hexa-thio triphenylene and hexa-thio azatriphenylene in two separate modelling/simulation steps. In the first we derived and parameterized a new GB molecular model for triphenylene-based compounds, by fitting the dispersive energies and the pair distributions calculated from the trajectories of an atomistic molecular dynamics (MD) simulation of chargeless hexa-octyl-thio triph-enylene (8HT-T). Then, by employing extensive Monte Carlo (MC) simulations we obtained the pressure-temperature phase diagram of the model, which shows crystalline, columnar, nematic, isotropic and gas regions. In the second step, we calculated, via *ab initio* methods, the electrostatic potential of hexa-thio triphenylene (HT-T) and hexa-thio azatriphenylene (HT-AT) cores and approximated them with two minimal sets of point charges [23]. We have studied, at a selected pressure, the variations induced on the phase behaviour of the two molecules by adding to the GB potential the point charges stemming from the HT-T and HT-AT electrostatic potential respectively.

## Results and Discussion

### Modelling

We modelled the total intermolecular potential as a sum of two additive terms: the first represents the dispersion forces and repulsive interactions between both the cores and the thioalkyl chains, and is intended to be a general model for a class of hexa-alkyl triphenylene derivatives with chain length ranging approximately from 5 to 9 atoms; the second is an electrostatic contribution, which is aimed at reproducing the Coulomb interactions between the molecular cores and which specifically depends on their chemical nature, here investigated for the triphenylene and azatriphenylene cases (Figure 1).

### Derivation of Gay-Berne parameters

The shape and the dispersive interactions of hexa-alkyl triphenylenes were described with an uniaxial repulsive-attractive GB [13] disc. The model considers only in an effective way the influence of molecular flexibility, by averaging over all the chain conformations.

*A*and

*B*as a function of the intermolecular distance vector

**r**

_{ AB }and of the orientation of the molecular orientation vectors

**u**

_{ A }and

**u**

_{ B }, with the help of six parameters: the exponents

*μ*and

*ν*which allow a tuning of the orientational dependence of the energy; the contact distances

*σ*

_{ ff }(face-face) and

*σ*

_{ ee }(edge-edge), in unit

*σ*

_{0}; the strength parameters

*ε*

_{ ff }and

*ε*

_{ ee }, in unit

*ε*

_{0}. The parameters

*σ*

_{ ff }and

*σ*

_{ ee }are respectively the lower and upper limit of the contact distance function

*σ*(

**u**

_{ A },

**u**

_{ B },

**r**

_{ AB }); the parameters

*μ*,

*ν*,

*ε*

_{ ff }and

*ε*

_{ ee }appear in the definition of the interaction function

*ε*(

**u**

_{ A },

**u**

_{ B },

**r**

_{ AB }) [13].

In most GB simulations these parameters are simply assumed a priori, or typically varied to look for generic trends in the phase diagram [12]. Here we wish instead to use parameters which are consistent with a coarse grained approximation to specific molecular organisations obtained from an atomistic MD simulation of a small sample of 8HT-T molecules in columnar phases (see Experimental Section).

using the atomic coordinates stored from the MD simulation at intervals of 10 picoseconds. We computed the intermolecular vectors **r**_{
AB
}of all pairs of molecules joining the centers of mass coordinates, and determined the two molecular orientations **u**_{
A
}and **u**_{
B
}as the eigenvectors of the molecular inertial frames corresponding to the maximum eigenvalues and calculated the corresponding GB energies. In order to obtain the best set of parameters for the system, we fitted all MD intermolecular potential energy points {**r**_{
AB
}, **u**_{
A
}, **u**_{
B
}, U_{
LJ
}} with the Gay-Berne expression by minimizing the total squared residuals ∑_{A,B;A≠B}[*U*_{
LJ
}(*A, B*) - *U*_{
GB
}(*A, B*)]^{2}. In the fitting procedure, the maximum GB energy was set to 10 kcal/mol and we also enforced the constraint of *ε*_{
ff
}> *ε*_{
ee
}. These two shape parameters were optimized using the Nelder-Mead symplex algorithm, while the *μ*, *ν*, *σ*_{
ff
}, *σ*_{
ee
}values were heuristically estimated by performing different series of minimization runs. The final parameters are: *μ* = 1, *ν* = 0, *σ*_{
ee
}= 19.25 Å, *σ*_{
ff
}= 3.75 Å, *ε*_{
ee
}= 9 kcal/mol, and *ε*_{
ff
}= 60 kcal/mol.

### Phase behaviour of charge-less model for triphenylene

Once the GB potential was parametrized, we studied by computer simulations the phase diagram of this coarse-grained molecular model (without charges: system Q0), by performing both MC cooling and heating scans.

*P** > 2, where

*P** ≡ $P{\sigma}_{0}^{3}/{\epsilon}_{0}$), the model exhibits a hexagonal columnar (Col

_{ h }) phase and at higher temperatures a nematic (Nem) phase. On the other hand, at low pressures and high temperatures, isotropic liquid (Iso) and gas phases are more stable. In the

*P** > 2 region the Nem-Col

_{ h }and Col

_{ h }-Cr transition temperatures increase as pressure increases (positive slope). Also the temperature ranges over which the nematic and the columnar phase are stable, widen with pressure. For

*P** < 0.5, the columnar phase disappears and the system shows a direct transition from a crystalline to an isotropic, liquid phase. The presence of the nematic phase is typical of rigid GB disks [17–21, 24], in which neither the Col

_{ h }-Iso transition can be favored by the flexibility of alkyl chains, nor the columnar phase can be stabilized, with respect to the nematic, by microsegregation of cores and chains. Indeed the heating phase diagram looks very similar to the ones reported by Caprion et al. [20] for GB discs with

*σ*

_{ ff }/

*σ*

_{ ee }= 0.2, and

*ε*

_{ ee }/

*ε*

_{ ff }= 0.1 and 0.2 (while here we have

*σ*

_{ ff }/

*σ*

_{ ee }= 0.1948, and

*ε*

_{ ee }/

*ε*

_{ ff }= 0.15). We identified the columnar phase as hexagonal with a spatial correlation length of about 15 molecules along the column axis (this was extimated from an additional MC simulation of a larger sample of 8000 molecules), and the crystalline one as interdigitated rectangular.

The sequence of cooling runs (blue dotted line in Figure 2B) shows a phase diagram similar to that obtained from heating runs, but with Cr-Iso and Cr-Col_{
h
}equilibrium curves shifted to higher *T** values (where *T** ≡ *kT*/*ε*_{0}) in the *P** < 2.5 region and to lower ones for *P** > 2.5. Also, the pressure stability range of Iso-Col_{
h
}co existence region appears to be strongly reduced with respect to the heating sequence. Indeed, at low pressure the transformation of the nematic phase into a Col_{
h
}structure did not proceed along a reversible path and showed significant hysteresis. Relevant hysteresis effects and even monotropic behavior have also been observed experimentally [4].

_{ h }transition temperature of 8H-TT (328 K) and compare it to the MC one (

*T** ≈ 0.5 at

*P** = 3.5) we obtain an estimate of the dimensioned value of

*ε*

_{0}: ${\epsilon}_{0}^{d}=k{T}_{\text{Cr}-{\text{Col}}_{\text{h}}}^{exp}/{T}_{\text{Cr}-{\text{Col}}_{\text{h}}}^{\ast}\simeq 1.3$ kcal/mol. To match the simulation and experiment temperature units, it is then sufficient to rescale it by a factor ${\epsilon}_{0}^{d}$/

*ε*

_{0}. We opted for this approach, taking for simplicity ${\epsilon}_{0}^{d}$ = 1 kcal/mol, having also considered that the experimental range of phase transition enthalpies (few kcal/mol in DSC measurements [4]) is well reproduced with this energy scale (see as example the enthalpy values in Figure 4D). Hence, using as unit dimensions

*ε*

_{0}= 1 kcal/mol and

*σ*

_{0}= 10 Å considering our definition for

*T** and

*P**, we obtained the following scaling factor:

*P*= 6.95

*P** MPa and

*T*= 500

*T** K, employed e.g. in Figure 2. Using these data, the limit simulation pressure

*P** = 2.5 turns out to be 17.5 MPa, and the intermediate

*P** = 3.5 turn out 24 MPa. While these values are greater than atmospheric pressure, they are not unrealistic, if we compare to the experimental pressures used in [4] for HT-T. Moving to the temperature scale and taking the intermediate pressure

*P** = 3.5 as reference, the Nem-Col

_{ h }transition temperature

*T** = 0.6 would correspond to a temperature of

*T*= 300 K.

### Derivation of electrostatic charges

^{+}geometry, by a direct fit of the molecular electrostatic potential. This procedure naturally determines partial charges for all atoms of the molecule, which could in turn be used to decorate the GB discs previously parametrized. However, exploring the whole phase diagram with such a large number of charges would be unfeasible, because of the computationally overwhelming number of electrostatic interacting pairs to be evaluated. Thus, following a charge-fitting methodology recently developed [23], we derived two minimal sets of point charges that reproduce the main features of both the molecular electrostatic potential and simultaneously the quadrupolar moment. The optimal number of charges suitably placed in the molecular frame resulted in being 12 for HT-T and 22 for HT-AT instead of the original 36 and 30 atomic partial charges (see Figure 5). We see at once from Figure 5 that a major difference of the two charge distributions is the presence of an effective central charge in the case of HT-AT while the core region is essentially neutral for HT-T. As in [23], the quality of the fit has been measured by comparing the intermolecular electrostatic energies of the full and reduced sets for a series of bimolecular configurations, with an overall root mean squared deviation of 0.26 and 0.06 kcal/mol respectively. Regarding the hypothesis that quadrupolar interactions are fundamental in determining the properties of the columnar phase [19], we believe that on the scale of short intermolecular distances, characteristic of the columnar arrangement (<4 Å), and more generally of condensed phases, the central multipolar description is not adequate and a local modelization of electrostatics is needed, justifying the choice of a point charges approach like the present one.

*σ*

_{0}, and with varying x and y components of their relative position vector. The twist angle between x molecular axes is fixed to 60°, which gives the maximum interaction [25]. The most striking difference is the presence of a deep minimum in the center of the GB+HT-T map, which, instead, is absent for the GB+HT-AT case, and this accounts for a reduced tendency of the molecules to stack on each other to form a columnar phase.

### Phase behaviour of charged discotics at *P** = 3.5

Since the phase diagram of system Q0 gives strong evidence of hysteresis betweeen heating and cooling runs, especially in the low pressure region, as observed experimentally [26], we concentrated our study on the effect of the different charge distributions at an intermediate pressure, *P** = 3.5, where the hysteresis effect is relatively negligible, and at temperatures corresponding to isotropic, nematic and columnar phases for system Q0. We therefore performed MC simulations of the GB model with the reduced point charges ${q}^{\ast}\equiv q/{({\epsilon}_{0}^{d}{\sigma}_{0})}^{1/2}$ representing the HT-T electrostatic potential (system Q1), and that for HT-AT (system Q2).

*T** dependence of average orientational order parameter ⟨

*P*

_{2}⟩, number density ⟨

*ρ**⟩ ≡ ${\sigma}_{0}^{3}$

*N*⟨1/

*V*⟩, adimensional energy ⟨

*U**⟩ ≡ ⟨

*U**⟩/

*ε*

_{0}and enthalpy ⟨

*H**⟩ ≡ ⟨

*H**⟩/

*ε*

_{0}per particle for systems Q0, Q1 and Q2. We first notice that all systems show an isotropic and a nematic phase, and that all average observables indicate a similar behavior of the Q1 and Q2 samples in these regions. Conversely the two systems with different charge distributions behave differently at low temperature. The impact of HT-T charges on the total energy is much more important than the one of HT-AT charges, inducing also a slight but significant increase in the orientational order and density at the lowest

*T** (see also Table 1). This is particularly evident at the Nem-Col

_{ h }transition (

*T** = 0.675), which for the HT-T case occurs at higher temperature than for system Q0, while it is absent in system Q2.

Observables at *P** = 3.5

system |
| ⟨ | ⟨ | ⟨ | ⟨ | phase |
---|---|---|---|---|---|---|

Q0 | 0.50 | 0.91 | 0.70 | -2.15 | - | Col |

0.55 | 0.89 | 0.67 | -2.00 | - | Col | |

0.60 | 0.81 | 0.62 | -1.50 | - | Nem | |

0.65 | 0.80 | 0.59 | -1.42 | - | Nem | |

0.70 | 0.78 | 0.58 | -1.34 | - | Nem | |

0.75 | 0.77 | 0.50 | -1.31 | - | Nem | |

0.80 | 0.69 | 0.55 | -1.26 | - | Nem | |

0.85 | 0.69 | 0.54 | -1.22 | - | Nem | |

0.90 | 0.68 | 0.52 | -1.18 | - | Nem | |

0.95 | 0.67 | 0.51 | -1.15 | - | Nem | |

1.00 | 0.59 | 0.49 | -1.14 | - | Nem | |

1.10 | 0.51 | 0.47 | -1.10 | - | Nem | |

1.20 | 0.4 | 0.45 | -1.00 | - | Nem | |

1.30 | 0.14 | 0.42 | -0.96 | - | Iso | |

1.40 | 0.07 | 0.40 | -0.90 | - | Iso | |

Q1 | 0.50 | 0.92 | 0.70 | -4.50 | -2.10 | Col |

0.55 | 0.90 | 0.68 | -4.18 | -1.98 | Col | |

0.60 | 0.88 | 0.66 | -3.92 | -1.78 | Col | |

0.65 | 0.87 | 0.63 | -3.83 | -1.75 | Col | |

0.70 | 0.79 | 0.59 | -2.54 | -1.00 | Nem | |

0.75 | 0.75 | 0.58 | -2.43 | -0.90 | Nem | |

0.80 | 0.69 | 0.55 | -2.29 | -0.90 | Nem | |

0.90 | 0.68 | 0.53 | -2.11 | -0.80 | Nem | |

1.00 | 0.62 | 0.49 | -1.70 | -0.55 | Nem | |

1.10 | 0.50 | 0.47 | -1.55 | -0.42 | Nem | |

1.20 | 0.40 | 0.45 | -1.39 | -0.38 | Nem | |

1.30 | 0.15 | 0.42 | -1.28 | -0.29 | Iso | |

1.40 | 0.07 | 0.40 | -1.20 | -0.27 | Iso | |

Q2 | 0.50 | 0.87 | 0.69 | -2.45 | -0.28 | Nem |

0.55 | 0.82 | 0.64 | -2.25 | -0.20 | Nem | |

0.60 | 0.76 | 0.61 | -1.66 | -0.14 | Nem | |

0.65 | 0.75 | 0.58 | -1.59 | -0.13 | Nem | |

0.70 | 0.74 | 0.57 | -1.50 | -0.12 | Nem | |

0.75 | 0.73 | 0.56 | -1.41 | -0.11 | Nem | |

0.80 | 0.68 | 0.54 | -1.36 | -0.10 | Nem | |

0.90 | 0.67 | 0.52 | -1.30 | -0.10 | Nem | |

1.00 | 0.60 | 0.49 | -1.19 | -0.09 | Nem | |

1.10 | 0.50 | 0.46 | -1.18 | -0.09 | Nem | |

1.20 | 0.39 | 0.45 | -1.10 | -0.09 | Nem | |

1.30 | 0.15 | 0.42 | -1.00 | -0.09 | Iso | |

1.40 | 0.08 | 0.40 | -0.90 | -0.08 | Iso |

*T** = 0.8) and isotropic (

*T** = 1.4) temperatures, while for the lowest temperature (

*T** = 0.6) of the HT-T discs (system Q1), columns are well defined and extend across the whole sample. On the contrary, even at this temperature system Q2 appears to be nematic. A more quantitative structural characterization of the phases is achieved by monitoring the positional order through the calculation of the radial distribution function g

_{0}(

*r*) and its second rank anisotropy ${g}_{2}^{+}(r)$. The former gives the probability of finding a molecule at a distance

*r*from the one at the origin, relative to the probability expected for a completely uniform distribution at the same density, while ${g}_{2}^{+}(r)$ =

*V*/(4

*πr*

^{2}

*N*)⟨

*δ*(

*r*-

*r*

_{ ij })

*P*

_{2}(cos

*β*

_{ ij })⟩

_{ ij }[27], where

*β*

_{ ij }is the angle between the intermolecular vector and the phase director, is a measure of the tendency of molecules to align in the face-face arrangement. In Figures 9A, 10A, we show the system Q1 and Q2 pair distribution functions for temperatures

*T** = 0.6, 0.8, 1.4. At the lowest temperature g

_{0}(

*r*) is very structured, indicating a high degree of positional order. The first three peaks at

*r*≈ 0.4, 0.8, 1.2

*σ*

_{0}correspond to neighbouring molecules stacked in the same column (face-face). The fourth peak is due to molecules (edge-edge) belonging to neighbouring columns (see also Figure 7). Less structured features are found in

*g*

_{0}(

*r*) at

*T** = 0.6 of system Q2, while at

*T** = 0.8, 1.4 the radial correlation functions show no positional ordering for both systems in study: edge-edge configurations are approximately as frequent as face-face. Consistently with the onset of the columnar phase the second rank anisotropy ${g}_{2}^{+}(r)$ (in Figures 9B, 10B) exhibits a rich structure even at large molecular separations: a close look at the curves shows three positive maxima for the intra-column neighbouring pairs (intermolecular vector parallel to the director) and a fourth negative peak corresponding to edge-edge molecules belonging to adjacent columns (intermolecular vector perpendicular to the director).

Also in the isotropic phase (*T** = 1.4) the HT-T system shows evidences of short range spatial and orientational correlation, i.e. also at high temperature we find couples of molecules in face-face arrangement. This configuration, like observed precedently for the contour plots, is favored by the HT-T charge distribution and strengthens the interaction between molecules at all temperatures (see *U** in Figure 4), with particular stabilization of the columnar phase.

## Conclusion

We derived a new parameterization of the GB potential by fitting the total LJ pair intermolecular potential from an atomistic MD simulation of hexa-alkyl-thio triphenylene. MC simulations show that the resulting molecular model exhibits isotropic, columnar and crystalline phases at low pressures, while at higher pressures (*P** > 2) a wide nematic phase appears and the columnar phase turns out enlarged. We believe this parametrization can be successfully employed as the basis for molecular level simulations of triphenylene molecules; notably different core sizes can be modelled by minor modifications of *σ*_{
ee
}and *σ*_{
ff
}parameters.

Adding the HT-T and HT-AT charges to the GB model at the selected pressure *P** = 3.5, we found that the intermolecular electrostatic potential among the cores is fundamental in stabilizing/destabilizing columnar phases; in particular the HT-T charge distribution stabilizes the columnar structure, while the HT-AT distribution suppresses its formation in favor of the nematic phase. This study therefore suggests that hexa-alkyl-thio azatriphenylenes do not form columnar phases due to the unfavorable charge distribution of the azatriphenylene core.

## Experimental

### MD simulation details

We ran atomistic MD simulations on systems of N = 40 8HT-T molecules, initially stacked in four columns of ten molecules each. We employed the AMBER-OPLS force field [28, 29] and we set all the point charges to zero, since the core-core electrostatic interaction have been taken into account in a separate specific model and the alkyl chains are approximately neutral.

On heating at ambient pressure, the 8HT-T solid compound melts from a crystalline phase to a columnar hexagonal phase at 328 K, and to the isotropic phase at 360 K [30]; in the simulation we used a temperature of 300 K and at a pressure of 3 atm to compensate the absence of electrostatic interactions and while remaining in a region of stability of the columnar arrangement. The simulations were run in the isothermal-isobaric ensemble with the Orac code [31]:, velocity scaling thermostat and isotropic Parrinello-Rahman barostat [32]: after an equilibration period of about 2 ns, we accumulated averages and trajectories for 2 ns more.

### MC simulation details

Systems Q0, Q1, Q2 were studied from MC cooling and heating scans at constant pressure, with temperature steps of Δ*T** ≤ 0.1. The equilibration stage of each simulation consisted of 300–400 thousand cycles and was followed by a production stage of 100 thousand cycles. Each MC cycle consists of trial displacements and reorientations of all the molecules and, with 10% probability, of a trial box side change. The simulation box was orthogonal with 3D periodic boundary conditions.

In systems with charges (systems Q1 and Q2) the long range electrostatic interactions were computed using the method proposed by Wolf *et al*. in [33], which is computationally more efficient than traditional Ewald summations. For the systems considered the compromise between accuracy and speed, with respect to standard Ewald method, is good provided that the truncation radius and the damping parameter *α* are adequately chosen (in this case *α* = 8). However, because of the high number of charges considered, simulation on systems Q1 and Q2 were run on pools of 8 processor using an MPI fortran77 program and employing a data replicated approach to distribute the energy computation between the threads.

In all experiments we considered a number of particles *N* = 1000, reduced pressures *P** ranging from 0.01 to 10, and a cutoff *r*_{
c
}of 3.5 *σ*_{0}, where the energy unit *ε*_{0} has initially been set to 30 kcal/mol and the distance unit *σ*_{0} to 10 Å. Phase transitions were located from the variation of the radial distribution function, of the second rank orientational order parameter ⟨*P*_{2}⟩ = $\u3008\frac{3}{2}{\left({z}_{i}\cdot Z\right)}^{2}-\frac{1}{2}\u3009$ and from the often small discontinuities in the enthalpy *H** = *U** + *P**/*ρ**.

## Declarations

### Acknowledgements

This work was supported by the EU Integrated Project NAIMO (NMP4-CT-2004-500355). We acknowledge prof. Yves H. Geerts of the Université Libre de Bruxelles for originally suggesting this investigation. We also thank Alan J. Miller for his package for Nelder Mead minimization [34].

## Authors’ Affiliations

## References

- Kumar S: Recent developments in the chemistry of triphenylene-based discotic liquid crystals. Liq Cryst. 2004, 31: 1037-10.1080/02678290410001724746.View ArticleGoogle Scholar
- Kumar S: Self-organization of disk-like molecules: chemical aspects. Chem Soc Rev. 2006, 35: 83-10.1039/b506619k.View ArticleGoogle Scholar
- Gramsbergen E, Hoving H, de Jeu W, Praefcke K, Kohne B: X-ray investigation of discotic mesophases of alkylthio substituted triphenylenes. Liq Cryst. 1986, 1: 397-View ArticleGoogle Scholar
- Maeda Y, Rao DSS, Prasad SK, Chandrasekar S, Kumar S: Phase behaviour of the discotic mesogen 2,3,6,7,10,11-hexahexyl thiotriphenylene (HHTT) under hydrostatic pressure. Liq Cryst. 2001, 28: 1679-10.1080/02678290110075084.View ArticleGoogle Scholar
- Roussel O, Kestemont G, Tant J, de Halleux V, Aspe RG, Levin J, Remacle A, Gearba IR, Ivanov D, Lehmann M, Geerts YH: Discotic liquid crystals as electron carrier materials. Molec Cryst Liq Cryst. 2003, 396: 35-10.1080/15421400390213177.View ArticleGoogle Scholar
- Palma M, Levin J, Lemaur V, Liscio A, Palermo V, Cornil J, Geerts Y, nad P, Samorì ML: Self-organization and nanoscale electronic properties of azatriphenylene based architectures: a scanning probe microscopy study. Advanced Materials. 2006, 18: 3313-10.1002/adma.200601437.View ArticleGoogle Scholar
- Gearba R, Lehmann M, Levin J, Ivanov D, Koch M, Barberá J, Debije M, Piris J, Geerts Y: Tailoring discotic mesophases: columnar order enforced with hydrogen bonds. Advanced Materials. 2003, 15: 1614-10.1002/adma.200305137.View ArticleGoogle Scholar
- Care CCM, Cleaver DJ: Computer simulation of liquid crystals. Reports on progress in physics. 2005, 11: 2665-10.1088/0034-4885/68/11/R04.View ArticleGoogle Scholar
- Wilson MR: Progress in computer simulations of liquid crystals. Int Rev Phys Chem. 2005, 24: 421-10.1080/01442350500361244.View ArticleGoogle Scholar
- Berardi R, Muccioli L, Zannoni C: Can nematic transitions be predicted by atomistic simulations? A computational study of the odd-even effect. ChemPhysChem. 2004, 5: 104-10.1002/cphc.200300908.View ArticleGoogle Scholar
- Muccioli L, Berardi R, Orlandi S, Ricci M, Zannoni C: Molecular properties and stacking of 1-substituted hexa-alkoxy-triphenylenes. Theor Chem Acc. 2007, [DOI 10.1007/s00214-006-0229-7]Google Scholar
- Zannoni C: Molecular design and computer simulations of novel mesophases. J Mat Chem. 2001, 11: 2637-10.1039/b103923g.View ArticleGoogle Scholar
- Gay JG, Berne BJ: Modification of the overlap potential to mimic a linear site-site potential. J Chem Phys. 1981, 74: 3316-10.1063/1.441483.View ArticleGoogle Scholar
- Berardi R, Orlandi S, Zannoni C: Antiphase structures in polar smectics liquid crystals and their molecular origin. Chem Phys Lett. 1996, 261: 357-10.1016/0009-2614(96)00983-9.View ArticleGoogle Scholar
- Berardi R, Ricci M, Zannoni C: Ferroelectric nematic and smectic liquid crystals from tapered molecules. ChemPhysChem. 2001, 2: 443-10.1002/1439-7641(20010716)2:7<443::AID-CPHC443>3.0.CO;2-J.View ArticleGoogle Scholar
- Emerson APJ, Luckhurst GR, Whatling SG: Computer-simulation studies of anisotropic systems XXIII. The Gay-Berne discogen. Molec Phys. 1994, 82: 113-10.1080/00268979400100094.View ArticleGoogle Scholar
- Bates MA, Luckhurst GR: Computer simulation studies of anisotropic systems XXVI. Monte Carlo investigations of a Gay-Berne discotic at constant pressure. J Chem Phys. 1996, 104: 6696-10.1063/1.471387.View ArticleGoogle Scholar
- Berardi R, Orlandi S, Zannoni C: Monte Carlo simulation of discotic Gay-Berne mesogens with axial dipole. J Chem Soc Faraday Trans. 1997, 93: 1493-1496. 10.1039/a607571a.View ArticleGoogle Scholar
- Bates MA, Luckhurst GR: Computer simulation studies of anisotropic systems XXIX. Quadrupolar Gay-Berne discs and chemically induced liquid crystal phases. Liq Cryst. 1998, 24: 229-10.1080/026782998207406.View ArticleGoogle Scholar
- Caprion D, Bellier-Castella L, Ryckaert JP: Influence of shape and energy anisotropies on the phase diagram of discotic molecules. Phys Rev E. 2003, 67 (4 Pt 1): 041703-10.1103/PhysRevE.67.041703.View ArticleGoogle Scholar
- Berardi R, Orlandi S, Zannoni C: Columnar and interdigitated structures from discotic mesogens with symmetric arrangment of dipoles. A Monte Carlo study. Liq Cryst. 2005, 32: 1427-10.1080/02678290500160696.View ArticleGoogle Scholar
- van de Craats AM: Charge transport in self-assembling discotic liquid crystalline materials. PhD thesis. 2000, Delft University of TechnologyGoogle Scholar
- Berardi R, Muccioli L, Orlandi S, Ricci M, Zannoni C: Mimicking electrostatic interactions with a set of effective charges. A genetic algorithm. Chem Phys Lett. 2004, 389: 373-10.1016/j.cplett.2004.03.119.View ArticleGoogle Scholar
- Berardi R, Cecchini M, Zannoni C: A Monte Carlo study of the chiral columnar organisations of dissymmetric discotic mesogens. J Chem Phys. 2003, 119: 9933-10.1063/1.1616913.View ArticleGoogle Scholar
- Lemaur V, da Silva Filho DA, Coropceanu V, Lehmann M, Geerts Y, Piris J, Debije MG, van de Craats AM, Senthilkumar K, Siebbeles LDA, Warman JM, Brédas JL, Cornil J: Charge transport properties in discotic liquid crystals: a quantum-chemical insight into structure-property relationships. J Am Chem Soc. 2004, 126: 3271-10.1021/ja0390956.View ArticleGoogle Scholar
- Rao DSS, Prasad SK, Bamezai RK, Kumar S: High pressure investigations on the phase behaviour of discotic liquid crystals. Molec Cryst Liq Cryst. 2003, 397: 143-View ArticleGoogle Scholar
- Berardi R, Emerson APJ, Zannoni C: Monte Carlo investigations of a Gay-Berne liquid crystal. J ChemSoc Faraday Trans. 1993, 89: 4069-10.1039/ft9938904069.View ArticleGoogle Scholar
- Cornell WD, Cieplak P, Bayly CI, Gould IR, Merz Jr KM, Ferguson DM, Spellmeyer DC, Fox T, Caldwell JW, Kollman PA: A second generation force field for the simulation of proteins and nucleic acids. J Am Chem Soc. 1995, 117: 5179-10.1021/ja00124a002.View ArticleGoogle Scholar
- Jorgensen WL, Maxwell DS, Tirado-Rives J: Development and testing of the OPLS all-atom force field on conformational energetics and properties of organic liquids. J Am Chem Soc. 1996, 118: 11225-11236. 10.1021/ja9621760.View ArticleGoogle Scholar
- Heiney PA, Fontes E, Jeu WHD, Riera A, Carroll P, Smith AB: Frustration and helicity in the ordered phases of a discotic compound. J Phys France III. 1989, 50: 461-View ArticleGoogle Scholar
- Procacci P, Paci E, Darden TA, Marchi M: ORAC: A molecular dynamics program to simulate complex molecular systems with realistic electrostatic interactions. J Comput Chem. 1997, 18: 1848-1862. 10.1002/(SICI)1096-987X(19971130)18:15<1848::AID-JCC2>3.0.CO;2-O.View ArticleGoogle Scholar
- Parrinello M, Rahman A: Polymorphic transitions in single crystals: A new molecular dynamics method. J App Phys. 1981, 52: 7182-7190. 10.1063/1.328693.View ArticleGoogle Scholar
- Wolf D, Keblinski P, Phillpot SR, Eggebrecht J: Exact method for the simulation of Coulombic systems by spherically truncated, pairwise r
^{-1}summation. J Phys Chem. 1999, 110: 8254-10.1063/1.478738.View ArticleGoogle Scholar - Software from Alan J. Miller. [http://users.bigpond.net.au/amiller]

## Copyright

This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.