Energetic and Entropic Contributions to the Landau–de Gennes Potential for Gay–Berne Models of Liquid Crystals

Abstract

The Landau–de Gennes theory provides a successful macroscopic description of nematics. Cornerstone of this theory is a phenomenological expression for the effective free energy as a function of the orientational order parameter. Here, we show how such a macroscopic Landau–de Gennes free energy can systematically be constructed for a microscopic model of liquid crystals formed by interacting mesogens. For the specific example of the Gay–Berne model, we obtain an enhanced free energy that reduces to the familiar Landau–de Gennes expression in the limit of weak ordering. By carefully separating energetic and entropic contributions to the free energy, our approach reconciles the two traditional views on the isotropic–nematic transition of Maier–Saupe and Onsager, attributing the driving mechanism to attractive interactions and entropic effects, respectively.

1. Introduction

There are two classical views on the isotropic-to-nematic (IN) transition in liquid crystals. One is based on the Maier–Saupe theory [1,2], where attractive interactions between the mesogens are supposed to be responsible for the alignment of the molecules. These attractive interactions are energetic in nature and therefore lead to the IN transition at sufficiently low temperatures, where they dominate over entropic effects. Liquid crystals where the transition is temperature-driven are called thermotropic. An alternative approach to the IN transition goes back to Onsager [3], where excluded volume interactions are identified as the driving force for orientational ordering. This athermal and purely entropic effect dominates for sufficiently elongated particles if the concentration is high enough. Since the IN transition is driven by density instead of temperature, those liquid crystals are termed lyotropic.

The Onsager approach is usually considered to be better founded since it has a well-defined limit of high-aspect ratio particles in the dilute regime. The corresponding density functional theories have been applied with good success to lyotropic liquid crystals [4,5]. The Maier–Saupe approach on the other hand lacks such a clear microsopic foundation. Nevertheless, due to its appealing simplicity it is used quite often in practice, especially for studying dynamical properties [6].

The Maier–Saupe and Onsager views emphasize different physical mechanisms, which typically are both present at the same time in real liquid crystals. Since both views are not mutually exclusive, there have been attempts to combine them [7]. In fact, the most macroscopic description of liquid crystal is provided by the Landau–de Gennes free energy as a function of the orientational order parameter [8]. As a macroscopic quantity, the Landau–de Gennes free energy is usually introduced phenomenologically as a low-order polynomial of the order parameter, in close analogy to the Ginzburg–Landau theory of phase transitions. Despite considerable advances in molecular simulations [9] and theoretical efforts [10,11], the long-standing question of how to derive the Landau–de Gennes free energy from a microscopic model has not been fully resolved yet. Since many length- and time-scales are present in typical liquid-crystalline systems, there is an urgent need for systematic coarse-graining approaches [12].

The Lebwohl–Lasher model [13,14] is one of the simplest models for liquid-crystals on the level of an interacting many-particle system including orientations. Excluded volume interactions are absent in this model since the molecule orientations are placed on a regular lattice. The model is constructed such that the Maier–Saupe theory is recovered in the mean-field limit. With the help of a systematic coarse-graining scheme, it is possible to obtain an enhanced Landau–de Gennes potential for the Lebwohl–Lasher model in non-polynomial form with physically admissible limits for strong orientational ordering [15].

Here, we study a more realistic model for liquid-crystals that includes translational as well as rotational degrees of freedom of the molecules. In particular, the interaction is chosen as the Gay–Berne potential [16], which is basically an appropriately scaled Lennard–Jones potential that allows two particles to approach each other closer parallel to the short axis compared with approaches along the long axis. The Gay–Berne potential is considered to be of a generic form and is therefore very widely used [17] for the simulation of liquid crystals in their isotropic, nematic and smectic phase. The phase diagram of the Gay–Berne model has been studied via computer simulations [18,19] and also by perturbation theory [20]. A suitable free energy functional has been proposed in the literature that allows to study the isotropic–nematic phase transition of the Gay–Berne model [5,21]. Nevertheless, the corresponding Landau–de Gennes potential still remains unknown. Here, we derive the Landau–de Gennes free energy for the Gay–Berne model following the general approach proposed in [15]. We use Monte-Carlo simulations in a generalized canonical ensemble and thermodynamic integration to work out separately the energetic and entropic contributions to the effective free energy.

2. From Micro to Macro

The point of departure is a microscopic system of N interacting particles, which is described by classical mechanics. We here consider a model of liquid crystals, where the individual molecules are described as rigid, anisotropic particles. The most prominent example for such systems is the Gay–Berne model [16] that we choose also here, but other models like [22] can be studied in the same manner. We want to emphasize that the approach followed here can also be applied to atomistic models of liquid crystalline molecules, for which the instantaneous orientation of the molecule or its rigid subunit is well-defined.

Let $\Gamma =\{{\mathbf{r}}_1,{\mathbf{u}}_1,\dots ,{\mathbf{r}}_N,{\mathbf{u}}_N\}$ denote the microstate of the system specified in terms of particle positions ${\mathbf{r}}_i$ and orientations ${\mathbf{u}}_i$. We restrict ourselves to particles that are rotationally symmetric. Therefore, the three-dimensional unit vectors ${\mathbf{u}}_i$ are sufficient to describe the orientation of the particles. In general, one needs to deal with the corresponding three Euler angles. The specific form of the Hamiltonian $H(\Gamma )$ completes the definition of the microscopic model system.

On a macroscopic level, we are interested only in the average orientation rather than those of all individual molecules. Therefore we specify the macrostate of the system by $(\rho ,T,\mathbf{Q})$, where ρ is the density, T the temperature, and Q the orientational order parameter. Due to head-tail symmetry of the molecules, Q cannot be a chosen as a vector quantity but as symmetric, traceless second-rank tensor [8]. For a given microscopic configuration Γ, the instantaneous value Π of the orientational order parameter tensor is given by

$$\mathbf{\Pi }(\Gamma )=\frac{1}{N}\sum _{j=1}^N({\mathbf{u}}_j{\mathbf{u}}_j-\frac{1}{3}\mathbf{1}),\mathbf{Q}=\langle \mathbf{\Pi }\rangle$$

(1)

where the brackets denote the thermal average with the probability density $\rho (\Gamma )$. The scalar Maier–Saupe orientational order parameter is defined as $S_2=(3/2)\mathbf{n}\mathbf{n}:\mathbf{Q}$, where n is the eigenvector corresponding to the largest eigenvalue of Q.

Systematic coarse graining from the microscopic to the macroscopic level of description is at the core of statistical thermodynamics. Significant progress in this field has been achieved in recent years (see e.g., [23,24]). We here follow the approach suggested in [15,25]. In agreement with the classical ideas of de Gennes, we make the crucial assumption that the order parameter is the only additional slow variable needed to describe the system on a macroscopic level [8]. With such a clear separation of scales, the relevant ensemble is given by the generalized canonical distribution ${\rho }_{\mathbf{\Lambda }}(\Gamma )\sim exp[-\beta H(\Gamma )-\mathbf{\Lambda }:\mathbf{\Pi }(\Gamma \left)\right]$, see Appendix A. The generalized canonical ensemble can be motivated by the maximum entropy principle (also termed “quasi-equilibrium approximation”), which has been used successfully for different kinetic models [26] including liquid crystals [27,28]. For Hamiltonian models of liquid crystals in particular, this ensemble has been employed in [15] in order to derive a generalized Landau–de Gennes potential for the Lebwohl–Lasher model. The Lagrange multiplier Λ ensures the prescribed average value $\mathbf{Q}={\langle \mathbf{\Pi }\rangle }_{\mathbf{\Lambda }}$, where the average is taken with ${\rho }_{\mathbf{\Lambda }}$. The normalization of the generalized canonical distribution is given by

$$\beta G(\mathbf{\Lambda })=-ln\int \mathrm{d}\Gamma exp[-\beta H(\Gamma )-\mathbf{\Lambda }:\mathbf{\Pi }(\Gamma \left)\right]$$

(2)

The potential G depends also on density and temperature. For ease of notation, however, we suppress the arguments ρ and T here and in the following. G can be interpreted as the generating function since the mean potential energy and orientational order parameter are given by

$$\begin{array}{ccc}\hfill E(\mathbf{\Lambda })& =& {\langle H\rangle }_{\mathbf{\Lambda }}=\frac{\partial \left(\beta G\right)}{\partial \beta }\hfill \\ \hfill \mathbf{Q}(\mathbf{\Lambda })& =& {\langle \mathbf{\Pi }\rangle }_{\mathbf{\Lambda }}=\frac{\partial \left(\beta G\right)}{\partial \mathbf{\Lambda }}\hfill \end{array}$$

(3)

The generalized free energy Ƒ is obtained as the Legendre transform of G with respect to Λ,

$$Ƒ\left(\mathbf{Q}\right)=G(\mathbf{\Lambda }\left(\mathbf{Q}\right))-k_{\mathrm{B}}T\mathbf{\Lambda }\left(\mathbf{Q}\right):\mathbf{Q}$$

(4)

The standard thermodynamic free energy is obtained as $F_0=Ƒ\left({\mathbf{Q}}_0\right)=G\left(0\right)$, where ${\mathbf{Q}}_0=\mathbf{Q}\left(0\right)$ is the equilibrium value of the order parameter for vanishing Λ. The generalized free energy Ƒ also serves as a thermodynamic potential, where the usual thermodynamic relationship $E=\partial \left(\beta Ƒ\right)/\partial \beta$ holds. In addition, the Lagrange multiplier obeys a similar relation as the dual variable to the orientational order parameter, $\mathbf{\Lambda }=-\partial \left(\beta Ƒ\right)/\partial \mathbf{Q}$, see Equation (21).

The generalized free energy can always be split as $Ƒ={Ƒ}^{\mathrm{id}}+{Ƒ}^{\mathrm{int}}$, where the contribution of the ideal orientational entropy ${Ƒ}^{\mathrm{id}}=-T{S}^{\mathrm{id}}$ in the absence of inter-particle interactions ($H=0$) has been worked out recently [29]. There, a rather accurate interpolation formula was found, $S^{\mathrm{id}}/k_{\mathrm{B}}=ln\left(\sqrt{Q}\right)+\frac{2}{3}(\sqrt{Q}-1)$, with $Q=1+9(I_3-I_2)$, with $I_2=(1/2)\mathrm{tr}\left({\mathbf{Q}}^2\right)$ and $I_3=\mathrm{tr}\left({\mathbf{Q}}^3\right)$ the second and third tensorial invariants of Q, respectively. Determining the remaining contributions of the interaction energy and excess entropy, ${Ƒ}^{\mathrm{int}}=\Delta E-T{S}^{\mathrm{exc}}$, is the main objective of this manuscript. In the following, we determine these contributions numerically from the free energy Ƒ. We calculate the latter via thermodynamic integration, which in the present case—thanks to the duality relation between Q and Λ—is simply given by

$$Ƒ\left(\mathbf{Q}\right)=F_0-k_{\mathrm{B}}T{\int }_{{\mathbf{Q}}_0}^{\mathbf{Q}}\mathbf{\Lambda }\left({\mathbf{Q}}^{\prime }\right):\mathrm{d}{\mathbf{Q}}^{\prime }$$

(5)

In the following, we will perform the thermodynamic integration in the isotropic state where ${\mathbf{Q}}_0=\mathbf{0}$.

3. Gay–Berne Model

In order to test the above ideas, we apply the approach to a specific model system of liquid crystals. It was shown in [15] that the approach outlined in Section 2 when applied to a simple lattice model of liquid crystals indeed leads to an enhanced Landau–de Gennes potential, which describes the numerical data rather accurately. Here, we apply this approach to a more realistic, off-lattice model system of ellipsoidal particles interacting via the Gay–Berne potential [16]. Except for kinetic energies, the Hamiltonian of the system is given by

$$H=\frac{1}{2}\sum _{i,j}{U}_{ij}({\mathbf{r}}_{ij},{\mathbf{u}}_i,{\mathbf{u}}_j)$$

(6)

where the pair potential depends not only on the distance vector ${\mathbf{r}}_{ij}$ but also on the relative orientation of the two particles,

$$U_{12}(\mathbf{r},{\mathbf{u}}_1,{\mathbf{u}}_2)=4\epsilon (\widehat{\mathbf{r}},{\mathbf{u}}_1,{\mathbf{u}}_2)[{\left(\frac{{\sigma }_0}{d(\mathbf{r},{\mathbf{u}}_1,{\mathbf{u}}_2)}\right)}^{12}-{\left(\frac{{\sigma }_0}{d(\mathbf{r},{\mathbf{u}}_1,{\mathbf{u}}_2}\right)}^6]$$

(7)

with $d(\mathbf{r},{\mathbf{u}}_1,{\mathbf{u}}_2)=r-\sigma (\mathbf{r},{\mathbf{u}}_1,{\mathbf{u}}_2)+{\sigma }_0$, where r is the distance between the centers of molecules 1 and 2. The orientation-dependent size parameter

$$\sigma (\mathbf{r},{\mathbf{u}}_1,{\mathbf{u}}_2)={\sigma }_0{\{1-\frac{\tilde{\chi }}{2}[\frac{{(\widehat{\mathbf{r}}·{\mathbf{u}}_1+\widehat{\mathbf{r}}·{\mathbf{u}}_2)}^2}{1+\tilde{\chi }{\mathbf{u}}_1·{\mathbf{u}}_2}+\frac{{(\widehat{\mathbf{r}}·{\mathbf{u}}_1-\widehat{\mathbf{r}}·{\mathbf{u}}_2)}^2}{1-\tilde{\chi }{\mathbf{u}}_1·{\mathbf{u}}_2}]\}}^{-1/2}$$

(8)

depends on the shape of the particles via the parameter $\tilde{\chi }=({\kappa }^2-1)/({\kappa }^2+1)$, where $\kappa ={\sigma }_{\mathrm{e}}/{\sigma }_{\mathrm{s}}$ is the length-to-breadth ratio of the particles. The energy term is usually written as

$$\epsilon (\widehat{\mathbf{r}},{\mathbf{u}}_1,{\mathbf{u}}_2)={\epsilon }_0{\epsilon }_1^{\mu }(\widehat{\mathbf{r}},{\mathbf{u}}_1,{\mathbf{u}}_2){\epsilon }_2^{\nu }({\mathbf{u}}_1,{\mathbf{u}}_2)$$

(9)

where

$$\begin{array}{ccc}\hfill {\epsilon }_2({\mathbf{u}}_1,{\mathbf{u}}_2)& =& {[1-{\tilde{\chi }}^2{({\mathbf{u}}_i·{\mathbf{u}}_j)}^2]}^{-1/2}\hfill \end{array}$$

(10)

$$\begin{array}{ccc}\hfill {\epsilon }_1(\widehat{\mathbf{r}},{\mathbf{u}}_1,{\mathbf{u}}_2)& =& 1-\frac{{\tilde{\chi }}^{\prime }}{2}[\frac{{(\widehat{\mathbf{r}}·{\mathbf{u}}_1+\widehat{\mathbf{r}}·{\mathbf{u}}_2)}^2}{1+{\tilde{\chi }}^{\prime }{\mathbf{u}}_1·{\mathbf{u}}_2}+\frac{{(\widehat{\mathbf{r}}·{\mathbf{u}}_1-\widehat{\mathbf{r}}·{\mathbf{u}}_2)}^2}{1-{\tilde{\chi }}^{\prime }{\mathbf{u}}_1·{\mathbf{u}}_2}]\hfill \end{array}$$

(11)

The parameter ${\tilde{\chi }}^{\prime }=({\kappa }^{\prime 1/\mu }-1)/({\kappa }^{\prime 1/\mu }+1)$ determines the anisotropy of the attractive forces, ${\kappa }^{\prime }={ϵ}_{\mathrm{s}}/{ϵ}_{\mathrm{e}}$ is the ratio of the energies for side-by-side over end-to-end configuration. The parameters ${\sigma }_0$ and ${ϵ}_0$ are the length and energy scale of the underlying Lennard–Jones potential for spherical particles, $\kappa ={\kappa }^{\prime }=1$.

We have chosen two sets of parameter values for which the Gay–Berne model has already been studied in the literature: $\kappa =3,{\kappa }^{\prime }=5,\mu =1,\nu =3$ in [30] and $\kappa =3,{\kappa }^{\prime }=5,\mu =2,\nu =1$ in [18,31]. In [30], the temperature-dependence of the model was investigated. There, the IN transition was found to occur at a reduced temperature $T^{*}=k_{\mathrm{B}}T/{ϵ}_0$ of $T^{*}\approx 3.6$ when the reduced density is fixed to ${\rho }^{*}=\rho {\sigma }_0^3=0.30$, see Figure 1. The concentration instead of the temperature was varied in [31]. There, the IN transition was located around ${\rho }^{*}\approx 0.32$ for a fixed temperature of $T^{*}=1.25$. The isotropic and nematic densities at coexistence for several temperatures were determined in [18]. Since the density changes at the isotropic–nematic transition are usually quite small, it is admissible to consider a fixed density. At even lower temperatures or higher densities, the Gay–Berne model shows a smectic phase, which, however, is beyond the scope of the present study. For clarity of presentation, in the following we focus on the results for the model parameters of [30].

We have performed Metropolis Monte-Carlo simulations of the Gay–Berne model (7) in the generalized canonical ensemble for systems containing $N=1000$ ellipsoidal particles. We have verified the simulation code by comparing various thermodynamic and structural quantities for the equilibrium case ($\mathbf{\Lambda }=0$) to literature data [30,31] (see also Figure 1).

For several sets of chosen Lagrange multipliers Λ, we evaluate the order parameter tensor $\mathbf{Q}(\mathbf{\Lambda })$ and other thermodynamic quantities of interest. The Lagrange multipliers act as orienting fields that induce a non-zero mean orientation even in the isotropic phase (3). In Figure 2, we show the orientational order parameter $Q_{11}$ for uniaxial Lagrange multipliers $\mathbf{\Lambda }=\lambda (\mathbf{n}\mathbf{n}-(1/3\left)\mathbf{1}\right)$ with n a unit vector. The orientational order increases with increasing ${\Lambda }_{11}=(2/3)\lambda$ until it saturates for $\lambda \gtrsim 10$. For weak fields and in the isotropic regime, a linear response relation holds, $\mathbf{Q}=\chi \mathbf{\Lambda }$ with the susceptibility χ. In the high-temperature, low-density limit, χ approaches the value for non-interacting rotors ${\chi }_0=-2/15$. For the present choice of parameters with a relatively high density, we do not observe this limit within the temperature range studied. Upon lowering the temperature or increasing the density, $\left|\chi \right|$ increases strongly near the isotropic–nematic transition.

Figure 1. Main frame: The Maier–Saupe orientational order parameter $S_2$ as a function of reduced temperature $T^{*}$ for the Gay–Berne model [30] from Monte-Carlo simulations of a system of $N=1000$ particles. The isotropic–nematic transition is clearly visible. The inset shows the mean potential energy per particle, $E/N$, also as a function of reduced temperature.

Figure 1. Main frame: The Maier–Saupe orientational order parameter $S_2$ as a function of reduced temperature $T^{*}$ for the Gay–Berne model [30] from Monte-Carlo simulations of a system of $N=1000$ particles. The isotropic–nematic transition is clearly visible. The inset shows the mean potential energy per particle, $E/N$, also as a function of reduced temperature.

Figure 2. (Color online) The main frame shows the order parameter $Q_{11}$ as a function of ${\Lambda }_{11}$ for uniaxial forms of the Lagrange multiplier. The model parameters of [30] have been chosen with different temperatures in the isotropic regime. The inset shows the initial susceptibility χ as a function of reduced temperature, corresponding to the linear relation shown as broken lines in the main frame.

Figure 2. (Color online) The main frame shows the order parameter $Q_{11}$ as a function of ${\Lambda }_{11}$ for uniaxial forms of the Lagrange multiplier. The model parameters of [30] have been chosen with different temperatures in the isotropic regime. The inset shows the initial susceptibility χ as a function of reduced temperature, corresponding to the linear relation shown as broken lines in the main frame.

We use not only uniaxial but also the general biaxial form for the Lagrange multipliers. In particular, we use the parametrization $\mathbf{\Lambda }=(\lambda +{\lambda }^{\prime }/2)(\mathbf{n}\mathbf{n}-(1/3)\mathbf{1})+{\lambda }^{\prime }(\mathbf{m}\mathbf{m}-(1/3)\mathbf{1})$, where n and m are two orthogonal unit vectors. We varied the strength of the biaxiality parameter ${\lambda }^{\prime }$ with respect to the unixial one λ within the range ${\lambda }^{\prime }/\lambda =0\dots 0.7$. Of special interest to us is the mean internal energy $E={\langle H\rangle }_{\mathbf{\Lambda }}$. In Figure 3, we plot the mean internal energy per molecule $E/N{k}_{\mathrm{B}}T$ as a function of the second invariant $I_2=(1/2)\mathrm{tr}\left({\mathbf{Q}}^2\right)$ for several choices of the Lagrange multipliers. The first result of this paper is the interesting finding that the mean internal energy can be described as a function of $I_2$ only, without explicit dependence on the third invariant $I_3$ even in the case of biaxial Lagrange multipliers. More quantitatively, we find that the mean internal energy is rather well described by the functional form

$$E(\rho ,T;\mathbf{Q})=E_0(\rho ,T)-u_1(\rho ,T)N{k}_{\mathrm{B}}T{I}_2-u_2(\rho ,T)N{k}_{\mathrm{B}}T{I}_2^2$$

(12)

The mean internal energy in the isotropic state $E_0$ is decreasing with decreasing temperature, showing a sudden drop at the isotropic–nematic transition (see inset in Figure 1). The second term in Equation (12) is the classical Maier–Saupe expression [1,2] for the mean interaction energy. The last term is a correction for stronger orientational ordering. A similar correction term has been proposed already in [32] in order to improve the comparison of the Maier–Saupe theory with experimental data. In agreement with the original ideas of Maier and Saupe, we observe that $u_1$ is positive and therefore favors parallel alignment of the molecules. For the model parameters suggested in [30], at a constant density ${\rho }^{*}=0.3$, the corresponding coefficients are given in Table 1. Since deviation of $E\left(I_2\right)$ from a linear form occurs only at rather strong orientational ordering, it is indeed justified to consider the quadratic term in Equation (12) as a correction, even though the magnitude of $u_1$ and $u_2$ are comparable, see Table 1. Moreover, we observe that $u_1$ and $u_2$ depend on the chosen model parameters, with $u_1$ increasing with decreasing temperature.

Figure 3. (Color online) Reduced mean internal energy per particle for the Gay–Berne model of [30] as a function of the second tensorial invariant of Q. From top to bottom, temperature decreases as $T=5.0,4.8,4.6,4.4,4.2$. Symbols represent the numerical values from Monte-Carlo simulations in the generalized canonical ensemble, while lines represent the fit to Equation (12).

Figure 3. (Color online) Reduced mean internal energy per particle for the Gay–Berne model of [30] as a function of the second tensorial invariant of Q. From top to bottom, temperature decreases as $T=5.0,4.8,4.6,4.4,4.2$. Symbols represent the numerical values from Monte-Carlo simulations in the generalized canonical ensemble, while lines represent the fit to Equation (12).

Table 1. Coefficients appearing in the Landau–de Gennes free energy Equation (14) and the mean internal energy Equation (12) for the model parameters used in [30]. Numbers in brackets denote the uncertainty in the last digit.

${\mathbf{\rho }}^{*}$	${\mathbf{T}}^{*}$	${\mathbf{s}}_{\mathbf{1}}$	${\mathbf{u}}_{\mathbf{1}}$	${\mathbf{u}}_{\mathbf{2}}$	${\mathbf{E}}_{\mathbf{0}}/\mathbf{N}{\mathbf{ϵ}}_{\mathbf{0}}$
0.3	5.0	1.20(1)	3.18(6)	2.0(3)	-0.167(2)
0.3	4.8	1.25(2)	3.24(8)	2.2(3)	-0.220(3)
0.3	4.6	1.21(1)	3.43(8)	2.0(3)	-0.275(2)
0.3	4.4	1.32(1)	3.44(8)	2.5(3)	-0.342(2)
0.3	4.2	1.26(1)	3.72(9)	1.7(4)	-0.420(3)
0.3	4.0	1.34(5)	3.73(8)	2.3(3)	-0.505(2)

Table 1. Coefficients appearing in the Landau–de Gennes free energy Equation (14) and the mean internal energy Equation (12) for the model parameters used in [30]. Numbers in brackets denote the uncertainty in the last digit.

${\mathbf{\rho }}^{*}$	${\mathbf{T}}^{*}$	${\mathbf{s}}_{\mathbf{1}}$	${\mathbf{u}}_{\mathbf{1}}$	${\mathbf{u}}_{\mathbf{2}}$	${\mathbf{E}}_{\mathbf{0}}/\mathbf{N}{\mathbf{ϵ}}_{\mathbf{0}}$
0.3	5.0	1.20(1)	3.18(6)	2.0(3)	-0.167(2)
0.3	4.8	1.25(2)	3.24(8)	2.2(3)	-0.220(3)
0.3	4.6	1.21(1)	3.43(8)	2.0(3)	-0.275(2)
0.3	4.4	1.32(1)	3.44(8)	2.5(3)	-0.342(2)
0.3	4.2	1.26(1)	3.72(9)	1.7(4)	-0.420(3)
0.3	4.0	1.34(5)	3.73(8)	2.3(3)	-0.505(2)

4. Landau–de Gennes Free Energy for Gay–Berne Model

Next, we evaluate the change in the generalized free energy from its isotropic value, $\Delta Ƒ=Ƒ-F_0$, by thermodynamic integration via Equation (5). The result is shown in Figure 4 for various temperatures. For the temperatures considered in Figure 4, the equilibrium corresponding to the minimum of $\Delta Ƒ$ is obtained for $I_2=0$, i.e., the isotropic state, in agreement with previous results (Figure 1). As expected, the minimum becomes shallower upon lowering the temperature towards the isotropic–nematic transition.

Figure 4. (Color online) Generalized free energy change per particle, $\Delta Ƒ/N{k}_{\mathrm{B}}T$, for the Gay–Berne model of [30] as a function of the second tensorial invariant of Q. From top to bottom, temperature decreases as $T=5.0,4.8,4.6,4.4,4.2,4.0,3.8$. The dashed vertical line represents the upper bound for the second invariant $I_2\le 1/3$.

Figure 4. (Color online) Generalized free energy change per particle, $\Delta Ƒ/N{k}_{\mathrm{B}}T$, for the Gay–Berne model of [30] as a function of the second tensorial invariant of Q. From top to bottom, temperature decreases as $T=5.0,4.8,4.6,4.4,4.2,4.0,3.8$. The dashed vertical line represents the upper bound for the second invariant $I_2\le 1/3$.

Taking advantage of the known form of the ideal orienational free energy for non-interacting particles ${Ƒ}^{\mathrm{id}}$ from [29], we subtract this term from $\Delta Ƒ$ and are left with the interaction contribution ${Ƒ}^{\mathrm{int}}=\Delta Ƒ-{Ƒ}^{\mathrm{id}}$. Using Equation (12) for the mean interaction energy, we can then determine the excess entropy $S^{\mathrm{exc}}$ from ${Ƒ}^{\mathrm{int}}=\Delta E-T{S}^{\mathrm{exc}}$, $\Delta E=E-E_0$, which is the only remaining unknown in the generalized free energy. Figure 5 shows the dimensionless excess entropy per particle, $S^{\mathrm{exc}}/N{k}_{\mathrm{B}}$, as a function of the second tensorial invariant $I_2$. It is interesting and surprising to note that, to a first approximation, $S^{\mathrm{exc}}$ can be well described by a simple function of $I_2$,

$$S^{\mathrm{exc}}(\rho ,T;\mathbf{Q})=s_1(\rho ,T)N{k}_{\mathrm{B}}\phi \left(I_2\right)$$

(13)

The slight deviations from a linear function can be rationalized by a mean-field approximation to Onsager’s excluded volume functional, $\phi \left(x\right)=1-\sqrt{1-3x}$ (see Appendix B and [33,34] for a derivation of $\phi$). The values of the coefficient $s_1$ obtained from a linear regression of the simulation data are shown in Table 1. As apparent from Figure 5, $s_1$ is only weakly increasing with decreasing temperature for the model parameters of [30]. Figure 5 demonstrates the entropic favoring of parallel alignment due to steric hindrance, i.e., even though orientational entropy is lost due to the alignment, translational entropy is increased since the available volume is larger for oriented particles [8]. We note that the simple form Equation (13) breaks down for very strong ordering. Since this regime is dominated by the divergence of the ideal orientational entropy $S^{\mathrm{id}}$, we expect deviations from Equation (13) to have a minor effect only.

Figure 5. (Color online) Dimensionless excess orientational entropy per particle, $S^{\mathrm{exc}}/N{k}_{\mathrm{B}}$, for the Gay–Berne model of [30] as a function of the second tensorial invariant of Q. The same color coding as in Figure 4 is used.

Figure 5. (Color online) Dimensionless excess orientational entropy per particle, $S^{\mathrm{exc}}/N{k}_{\mathrm{B}}$, for the Gay–Berne model of [30] as a function of the second tensorial invariant of Q. The same color coding as in Figure 4 is used.

Putting the results from the simulations of the Gay–Berne model together, we arrive at the Landau–de Gennes potential,

$$\begin{array}{ccc}\hfill Ƒ(\rho ,T,\mathbf{Q})& =& F_0(\rho ,T)+{Ƒ}^{\mathrm{id}}(T,\mathbf{Q})+{Ƒ}^{\mathrm{int}}(\rho ,T,\mathbf{Q})\hfill \\ \hfill {Ƒ}^{\mathrm{id}}(T,\mathbf{Q})& =& -T{S}^{\mathrm{id}}\left(\mathbf{Q}\right)=-N{k}_{\mathrm{B}}T[ln\left(\sqrt{Q}\right)+\frac{2}{3}(\sqrt{Q}-1)]\hfill \\ \hfill {Ƒ}^{\mathrm{int}}(\rho ,T,\mathbf{Q})& =& -N{k}_{\mathrm{B}}T[u_1(\rho ,T)I_2+u_2(\rho ,T)I_2^2+s_1(\rho ,T)\phi \left(I_2\right)]\hfill \end{array}$$

(14)

where Q was defined above. Equation (14) is the main result of the present study, which provides the most economical macroscopic description of the Gay–Berne model system.

5. Isotropic–Nematic Transition from Landau–de Gennes Free Energy

The standard Landau–de Gennes form of the effective free energy is obtained by expanding Ƒ from Equations (14) for weak orientational ordering,

$$Ƒ/N{k}_{\mathrm{B}}T\approx F_0/N{k}_{\mathrm{B}}T+\frac{1}{2}(\frac{15}{2}-c_1)\mathrm{tr}\left({\mathbf{Q}}^2\right)-\frac{75}{14}\mathrm{tr}\left({\mathbf{Q}}^3\right)+\frac{1}{4}(\frac{3825}{196}-c_2){\left[\mathrm{tr}\left({\mathbf{Q}}^2\right)\right]}^2+\mathcal{O}\left({\mathbf{Q}}^5\right)$$

(15)

where $c_1=u_1+(3/2)s_1$ and $c_2=u_2+(9/2)s_1$. From Equation (15) we recognize the familiar scenario of a first order phase transition within Landau theory [8]. Different from the standard expression, the coefficient of the fourth order term is slightly density- and temperature-dependent, due to the deviation of the mean internal energy from the Maier–Saupe expression for strong orientational ordering and nonlinear entropic effects. In order to describe stable biaxial phases, a sixth order expansion of the free energy was found to be necessary [35], which clearly is contained within our approach. Figure 6 shows the full Landau–de Gennes free energy Ƒ, Equation (14), as a function of the Maier–Saupe orientational order parameter $S_2$ for the model parameters of [30]. While the overall behavior is characteristic for a first order transition, clear deviations from a polynomial behavior are seen due to the singular behavior of Ƒ near the full oriented states at $S_2=-1/2$ and $S_2=1$. Within the isotropic regime, the analytical formula agrees nicely with the simulation results from thermodynamic integration. Extrapolating the model parameters given in Table 1 to lower temperatures, we include in Figure 6 also the graphs for Ƒ at the isotropic–nematic transition, which we estimate as $T_{\mathrm{IN}}\approx 3.52$, and at an even lower temperature in the nematic phase. For the value of $S_2$ at the transition we find $S_2\left(T_{\mathrm{IN}}\right)\approx 0.59$, in good agreement with previous results (see Figure 1).

Figure 6. (Color online) The dimensionless Landau–de Gennes free energy per particle, $Ƒ/N{k}_{\mathrm{B}}T$ from Equation (14), as a function of the Maier–Saupe orientational order parameter $S_2$, for the Gay–Berne model of [30]. Temperature decreases from top to bottom as $T^{*}=5.0,4.0,3.52,3.3$. Symbols denote the result of Monte-Carlo simulation in the generalized canonical ensemble in the isotropic state.

Figure 6. (Color online) The dimensionless Landau–de Gennes free energy per particle, $Ƒ/N{k}_{\mathrm{B}}T$ from Equation (14), as a function of the Maier–Saupe orientational order parameter $S_2$, for the Gay–Berne model of [30]. Temperature decreases from top to bottom as $T^{*}=5.0,4.0,3.52,3.3$. Symbols denote the result of Monte-Carlo simulation in the generalized canonical ensemble in the isotropic state.

In the present model, the isotropic–nematic transition occurs when the combined energetic and excess-entropic gain for orientational ordering overwhelms the loss of ideal orientational entropy. Therefore, Onsager’s excluded volume interaction and Maier–Saupe’s attraction both play a role on the macroscopic level of this model, their mutual importance depending on the chosen state point. Since we find that the dominant energetic contribution $u_1$ grows with decreasing temperature in a similar manner as the entropic contribution, we conclude that even in the case where the IN transition is temperature-driven, entropic effects can play an important role.

For the case of hard ellipsoids, it was found by extensive simulations that a minimum axis ratio is needed in order to allow for IN transition [36]. This finding is consistent with the present approach since there is no energetic contribution to the free energy for hard particles, so that entropic effects have to compensate for that. The density- and temperature-dependence of the coefficient $a(\rho ,T)=15/2-c_1(\rho ,T)$ of the quadratic term in Equation (15) is crucial for the isotropic–nematic phase transition. The condition $a=0$ defines the pseudocritical point where the isotropic state becomes unstable. From $\mathbf{\Lambda }=-\partial \left(\beta Ƒ\right)/\partial \mathbf{Q}$, Equation (21), we find the simple relation $\chi =1/a$ in the isotropic regime, with the susceptibility χ of orientational ordering in response to a weak applied field Λ, $\mathbf{Q}=\chi \mathbf{\Lambda }$. Therefore, the susceptibility χ diverges not at the IN transition but at a lower temperature, where the isotropic state loses its stability (see inset to Figure 2). Note that the susceptibility is in general anisotropic in the nematic regime and therefore has to be represented by a tensorial quantity.

6. Conclusions

We have demonstrated a systematic method in order to derive the macroscopic Landau–de Gennes free energy of nematics from an underlying microscopic model. The method relies on the generalized canonical distribution and corresponding thermodynamic integration. From Monte-Carlo simulations in the generalized canonical ensemble, we find explicit expressions for the energetic and entropic contributions to the effective free energy. The method extends earlier work along the same lines on the Lebwohl–Lasher model [15] to off-lattice models including excess entropy due to translational degrees of freedom. For the special case of ellipsoidal particles interacting via the Gay–Berne potential with the parameters of [30], we find that the excess entropy is an additional driving force to energetic attraction for the isotropic–nematic transition, where the latter dominates for this temperature-driven case.

Besides addressing such questions, the method presented here can be useful in large-scale studies of defect structures in liquid crystals or hybrid systems, such as colloidal particles dispersed in a liquid crystal [37]. The multi-scale nature of these systems prevents a detailed simulation of liquid-crystalline mesogens such that a coarse-grained modeling has to be employed. While so far the standard Landau–de Gennes free energy was used in these studies with little or no connection to the underlying microscopic model, our approach allows to use a more realistic macroscopic free energy that is determined by the properties of the microscopic model. Since also surface tension effects in liquid crystals are ruled by the effective free energy [38], a consistent treatment involving interfaces could benefit from the present study. For studying phase ordering kinetics, it was found experimentally that thermal diffusivity changes the scaling laws for domain growth in liquid crystals [39]. Working within a non-equilibrium thermodynamics framework, our approach should also be useful in addressing this issue. Finally, the Landau–de Gennes free energy is routinely used in phenomenological models for nematodynamics [8]. There, especially in strong flow conditions, a quartic form of the free energy is clearly insufficient to keep the order parameter within the physically admissible range [34] and an enhanced Landau–de Gennes potential is very desirable.