Fundamental measure theory for smectic phases: Scaling behavior and higher order terms René Wittmann, Matthieu Marechal, and Klaus Mecke Citation: The Journal of Chemical Physics 141, 064103 (2014); doi: 10.1063/1.4891326 View online: http://dx.doi.org/10.1063/1.4891326 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/141/6?ver=pdfcov Published by the AIP Publishing Articles you may be interested in Effect of polydispersity and soft interactions on the nematic versus smectic phase stability in platelet suspensions J. Chem. Phys. 134, 124904 (2011); 10.1063/1.3570964 Phase behavior of ionic liquid crystals J. Chem. Phys. 132, 184901 (2010); 10.1063/1.3417384 Smectic ordering in athermal systems of rodlike triblock copolymers J. Chem. Phys. 127, 154902 (2007); 10.1063/1.2787009 A novel orientation-dependent potential model for prolate mesogens J. Chem. Phys. 122, 024908 (2005); 10.1063/1.1830429 The phase behavior of a binary mixture of rodlike and disclike mesogens: Monte Carlo simulation, theory, and experiment J. Chem. Phys. 119, 5216 (2003); 10.1063/1.1598432

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

THE JOURNAL OF CHEMICAL PHYSICS 141, 064103 (2014)

Fundamental measure theory for smectic phases: Scaling behavior and higher order terms René Wittmann, Matthieu Marechal, and Klaus Mecke Institute for Theoretical Physics, Universität Erlangen-Nürnberg, Staudtstr. 7, 91058 Erlangen, Germany

(Received 3 June 2014; accepted 15 July 2014; published online 8 August 2014) The recent extension of Rosenfeld’s fundamental measure theory to anisotropic hard particles predicts nematic order of rod-like particles. Our analytic study of different aligned shapes provides new insights into the structure of this density functional, which is basically founded on experience with hard spheres. We combine scaling arguments with dimensional crossover and motivate a modified expression, which enables an appropriate description of smectic layering. We calculate the nematic– smectic-A transition of monodisperse hard spherocylinders with and without orientational degrees of freedom and present the equation of state and phase diagram including these two liquid crystalline phases in good agreement with simulations. We also find improved results related to the isotropic–nematic interface. We discuss the quality of empirical corrections and the convergence towards an exact second virial coefficient, including higher order terms. © 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4891326] I. INTRODUCTION

Density functional theory (DFT)1 is a powerful tool to study inhomogeneous liquids. One of the most successful functionals, Rosenfeld’s Fundamental measure theory (FMT),2 gives rise to the Percus-Yevick equation of state (EOS)3 for the fluid phase of hard spherical particles, as it is based on scaled particle theory.4 The dimensional crossover to extremely confined systems allowed a subsequent finetuning of the functional, which made the quantitative description of the solid state accessible.5–8 Starting with zerodimensional cavities, the functional can be systematically rederived such that the free energy in lower dimensions is recovered for sharply peaked density distributions.7 The final version of Tarazona8 is constructed to reproduce the PercusYevick direct correlation function, which is a basic feature of the original functional.2 A generalized dimensional crossover approach for anisotropic hard particles can be found in Ref. 9. Inhomogeneous phases of anisotropic particles were studied within various fundamental measure approaches. All earlier functionals either focus on some specific class of anisotropic particles or are restricted to discrete orientations: Functionals for parallel squares and cubes were constructed from the zero-dimensional limit.10–13 Parallel hard cylinders (HC)14, 15 can be described by demanding a correct crossover to two-dimensional FMT functionals.2, 7 The generalization of hard cubes to parallelepipeds16 also served as a basis for the study of rods and plates in the Zwanzig model.17, 18 The restriction to infinitely thin particles gave rise to functionals for mixtures of spheres and rods19, 20 and later also platelets.21 A more general functional for arbitrary particle orientations was derived from an interpolation between the Rosenfeld functional2 and the Onsager approach,22 which makes it only suitable for prolate particles.23 The Rosenfeld functional itself is unable to describe nematic phases

0021-9606/2014/141(6)/064103/19/$30.00

when applied to anisotropic bodies.24, 25 Extended deconvolution fundamental measure theory (edFMT), a generalization which involves additional, anisotropy-sensitive, geometrical expressions, has been successfully applied to the isotropic– nematic transition26, 27 and the inhomogeneous interface28 of hard spherocylinders (HSC). Other applications are polyhedra near a hard wall29 and parallel hard cubes with rounded edges.30 The phase diagram of HSC is well understood by means of computer simulation.31–33 In particular, the nematic– smectic-A transition is of first order even for infinitely long rods.33 Poniewierski used the virial expansion up to third order, to confirm this result theoretically.34 Onsager’s second virial approximation22 alone, is not sufficient to predict a first order transition.34 Simulations also show that the coexisting packing fractions decrease with increasing aspect ratio.32 This has been an issue for a long time due to contradictory DFT results. The wrong positive slope was found with the weighted density approach by Poniewierski and co-workers35–38 while the functional by Somoza and Tarazona39–42 yields the opposite behavior. A tricritical point is found within both approaches for considerably different aspect ratios. Several more simplistic extended Onsager theories based on the Parsons-Lee approach43, 44 were applied in Refs. 45 and 46. Graf and Löwen proposed a modified weighted density approximation47 and an approach based on a cell model.48 All functionals accordingly predict the location of the isotropic–nematic–smectic-A triple point for aspect ratios roughly between 2.7 < l < 4 in agreement with simulations.32 Further simplifications can be achieved when the orientational degrees of freedom are neglected: Perfectly aligned HSC are an artificial but easily accessible system. MonteCarlo simulations49–51 show that the nematic–smectic transition is of second order. In an early study, this transition was located by means of bifurcation analysis for the virial

141, 064103-1

© 2014 AIP Publishing LLC

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-2

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

expansion.52 As a first step, more advanced density functionals were applied to parallel HSC.36, 39, 40, 45 A common finding is the decreasing nematic–smectic transition density with increasing aspect ratio, while quantitative results differ, especially between simulation and DFT.15 Additional orientational degrees of freedom increase the complexity significantly. A common assumption in theory is the orientational decoupling approximation, which states that the orientational degrees of freedom do not depend on the spatial coordinates.35, 37, 41, 45, 47, 48 Generally, this is considered to be a crude approximation, as a strong spatial modulation of the orientational distribution function is observed in both simulations and Onsager DFT.53 The additional assumption that the smectic orientational distribution is equal to the nematic at the respective density is not sufficient to describe a first order transition.35, 37 A free order parameter within the same distribution, instead, turned out to be a successful with deviations smaller than 0.2% to a minimization without orientational decoupling.42 The reason is that the low-density, less ordered regime between the smectic layers does not contribute significantly to the free energy. In Sec. II we review different versions of edFMT and provide a generalized dimensional crossover study. We present an extended scaling analysis in Sec. III, in order to ensure a valid description of arbitrary shapes at higher densities. We motivate an analytic approach to perfectly aligned systems. The main part in Sec. IV is dedicated to the smectic phase of HSC. After an analytic study of perfect alignment, we consider the phase diagram of orientational disordered HSC including isotropic, nematic and smectic-A phases. We find the triple point for an aspect ratio between 2.9 < l < 2.95 and a second order transition for l > 13. We end up with a discussion in Sec. V regarding the consistency of current approximations. In Appendix A, we present an analytic approach to calculate the smectic free energy of infinitely long HSC, using a scaled functional based on results shown in Appendix B.

satisfactory quantitative description.28 The quite accurate results with a more ad hoc functional41, 42, 62 suggest that troubleshooting and improvements to edFMT should have a high priority. To this end, we will stick to monodisperse bulk systems of simple body shapes in this work. We will show that we can choose restricted systems, which do not require any empirical modifications. Within this ideal framework, we will test individual components of the functional and approximations, which we derive below. Given a density functional     β ([ρ], μ) = dr id (r) + ex (r) − βμ dρ(r,  )

II. DENSITY FUNCTIONAL FOR ANISOTROPIC HARD BODIES

for anisotropic hard bodies. The inhomogeneous number density ρ(r,  ) depends on position r and orientation  . The geometric weight functions   ˆ − |r| , ω(3) (R) = |R(R)|

First, we review the state of the art of the edFMT functional. In Sec. II A, we consider some variants to include the second virial coefficient with different accuracy and discuss the functional forms for higher densities, proposed for hard spheres. We re-derive the functional from a generalized dimensional crossover point of view in Sec. II B, which holds for an arbitrarily shaped body. A. Extended fundamental measure theory

Rosenfeld’s fundamental measure theory provides the most founded approach to the phase behavior of arbitrary mixtures of hard spherical bodies bounded by any external potential.2 The generalization to arbitrarily shaped convex bodies, edFMT, shares both, a straightforward derivation from first principles and the simple functional structure, where all building blocks only depend on the geometry of one single body.26, 27 It is logical to assume that edFMT has the potential to reach equivalent success. However, current approaches use approximations with empirical inputs and still lack of a

(1) of the grand potential, a straightforward minimization δ/δρ = 0 provides an adequate description of the equilibrium state.1 Assuming orientational decoupling, the inhomogeneous number density ρ(r,  ) = ρ(r)g( ) = ρh(r)g( )

(2)

can be factorized into a spatial modulation ρ(r) of the average number density ρ = ηv, described with a modulation function h(r) of position r, and a normalized distribution g( ) of orientations  . The packing fraction of bodies with volume v is denoted by η. The contribution of the non-interacting ideal gas to the free energy density reads    id (r) = ρ(r) ln (ηh(r)) − 1 + d g( ) ln g( ) , (3) where we chose 3 = v for the thermal wavelength . As described in Refs. 26 and 27 we can use the weighted densities    nν (r) = dr (4) d ρ(r  ,  ) ω(ν) (r − r  ,  ) , to write down the excess free energy density φ2 φ3 ex ({ni }) = −n0 ln(1 − n3 ) + + 1 − n3 (1 − n3 )2

ω(2) (R) =

(5)

ˆ − |r|) δ(|R(R)| , ˆ r n(R)ˆ

ˆ H (R) ω(2) (R), 4π ˆ K(R) ω(0) (R) = ω(2) (R), 4π − → ˆ ω(2) (R), ω (2) (R) = n(R) ω(1) (R) =

(6)

ˆ − H (R) → − → ω (2) (R), ω (1) (R) = 4π ˆ . . . n (R) ˆ ω(2) (R), ωα(2)...α (R) = n (R) 1

r

ωα(1)1 ...αr (R) =

α1

αr

ˆ κ(R) ˆ . . . n (R) ˆ ω(2) (R) nα (R) αr−2 1 4π   ˆ vαI (R) ˆ − vαI I (R) ˆ vαI I (R) ˆ × vαI r−1 (R) r r−1 r

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-3

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

characterize the shape of the bodies in terms of mean H = 12 (κ I + κ II ), Gaussian K = κ I κ II , and deviatoric κ = 12 (κ I − κ II ) curvatures, the outward surface normal n and the directions vI/II of principal curvatures κ I/II on the surface parametrized by R. All quantities depend on position and orientation, combined in R = (r,  ), the radial unit vector ˆ = (ˆr,  ). rˆ = r/|r|, and respectively, R There are different expressions for the second φ 2 and third term φ 3 which all need to reduce to the Rosenfeld functional2 when applied to the homogeneous hard sphere fluid. The second virial coefficient B2 ρ 2 = n0 n3 + φ2

(7)

can be described exactly through the second term → → φ2[∞] = n1 n2 − − n 1− n2−

∞   [r] [r] n1 n2

(8)

r=2

of the functional, which involves a series of contraction products 3 3    [r] [r] α ...α α ...α n1 n2 := ... n1 1 r n2 1 r α1 =1

1

2

1

2

26–28

for the description of anisotropic bodies. The missing higher order terms are corrected with the empirical ζ parameter which is sufficient for a qualitative characterization of the nematic phase of HSC. A semi-empirical effort, which assigns orientation dependence to a running parameter ζ (S), was made to improve upon a constant ζ correction.27 Note that non-spherical bodies with small anisotropy as, e.g., dumbbells9 are adequately described with ζ = 0, i.e., the original Rosenfeld proposal.24, 25 However, the rank two tensors [r ] in φ2 t vanish for the homogeneous fluid of bodies with cubic symmetry, even if they are aligned.30 Therefore, fourth order tensors are essential in this case. We propose a generalized correction (ζ [nt ])

φ2

→ → n2− = n1 n2 − − n 1−

φ3(TR) =

nt n  

 ζm[n] n[2m] n[2m] 1 2

φ3(T) =

of higher order terms, truncated after rank nt tensors. This expansion does not incorporate tensors of odd rank, as the

(14)

3 − → → → → → n 2− n 2 − n2 − n 2− n2 n T2 ← 16π  →3  →2  −Tr ← n + n Tr ← n 2

2

(15)

2

introduced by Tarazona restores the Percus-Yevick result8 This final expression is considered to be the most sophisticated, concerning the freezing of hard spheres. Therefore, it was introduced to edFMT as well.26, 27 For anisotropic bodies, we expect the phase behavior to critically depend on the expressions for φ 2 and φ 3 which contain tensorial weighted densities and thus crucial information about anisotropy. In comparison, a generalization of Eq. (5), which adjusts prefactors fn (n3 ) to n0 , φ 2 , and φ 3 to obtain the Carnahan-Starling EOS54 for hard spheres,55, 56 does not improve the isotropic EOS of HSC significantly.28 B. Dimensional crossover for arbitrary shapes

One way to derive a functional for FMT, is to carry out a crossover to zero dimensions6, 7 for a system in the grand canonical ensemble. First, one considers a cavity, which is so small that at maximum one particle fits in. Subsequently, the cavity is shrunk while increasing the chemical potential, such that the average number of particles N ≤ 1 is fixed. In the limit, where only one configuration is accessible to the particle, the functional can be obtained analytically by comparing to the exact free energy, which is known for any cavity for which N ≤ 1. As shown by Tarazona and Rosenfeld,7 the resulting functional is equal to the first term in Eq. (5). Considering two overlapping cavities such that N ≤ 1 results in the respective first two terms with the exact φ2[∞] from Eq. (8) and, finally, considering three cavities such that N ≤ 1 gives a φ 3 which has the form

(12)

n=1 m=1

 →3  3  3 → n − 3n2 Tr[← n 22 ] + 2Tr ← n 2 16π 2

by Tarazona and Rosenfeld.7 It has the right dimensional crossover properties but gives rise to a different direct correlation function for hard spheres. The subsequent version

(9)

to Rosenfeld’s original proposal25 to decompose the Mayer ffunction for the interaction of two bodies Bi and Bj , where ∂ denotes the boundary.26, 27 The expansion of the denominator 1 + ni nj allows the factorization into independent contributions of one body, namely, the tensorial weight functions [r ] in Eq. (6). However, the truncated versions φ2 t of Eq. (8) at maximal tensor rank rt converge slowly for HSC.27 So far → primarily rank two tensors ← n ν ≡ n[2] ν have been used in the truncated second term (ζ ) → → → → φ =n n −− n − n ← (11) n − ζ Tr[← n ] 1 2

(−1)n+m 4−n (4n + 1)(2n + 2m)! 2n(n + 1)(4n2 − 1)(n + m)!(n − m)!(2m − 2)! (13) in Appendix E which implies the first order correction (ζ ) ζ = ζ1[1] = 5/4 for φ2 in Eq. (11). We will stick to the use of (ζ ) φ2 in large parts of this work for the reason of strongly increasing complexity of other approaches for inhomogeneous problems. A more comprehensive study of the complete tensor series, Eq. (8), will be part of Secs. III B and IV B which are dedicated to parallel HSC. The original third term of Rosenfeld2 was substituted with the expression ζm[n] =

αr =1

of rank r tensorial weighted densities n[r] ν . This series is directly related to the additional term 2  2   ds κ i vIi nj − vII i nj κ Iij = (10) (1 + ni nj ) |ni × nj | ∂Bi ∩∂Bj 4π

2

coefficients go back to the mutual excluded volume of HSC, studied in Appendix D. Explicitly, we calculate

φ3 :=



3 

di dri ρ(ri , i )ω(2) (r − ri , i )



i=1

×K3 ({nj (r − rj , j )}j =1,2,3 ) .

(16)

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-4

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

The kernel, which depends only on the normal vectors {nj (r − ri , j )}j =1,2,3 of the three surfaces, reads7, 9 (I)



K3 ({nj }) =

 |n1 · (n2 × n3 )| 2π − ψi−1,i,i+1 24π i=1 3

(II)

(IIIa)

(IIIb)

(IV)

(17)

where  ψij k := arccos (ni × nj ) · (nj × nk ) ,

(18)

for i, j, k ∈ {1, 2, 3}, while ψ 012 := ψ 312 and ψ 234 := ψ 231 . To write this third term in terms of weighted densities,8 it has to be expanded in the normal vectors. The only expansion that (i) leads to a φ 3 with only scalar, vector and second order tensorial weight functions (ii) is isotropic under simultaneous rotation of the three bodies and (iii) vanishes when two cavities are considered is app

K3

=

ξA ξ

(1 − ni · nj ), (19) [n1 · (n2 × n3 )]2 + B 16π 16π {ij }

Note, that this expansion contains two arbitrary parameters, which do not a priori give rise to a Rosenfeld-like expression, which yields the Percus-Yevick EOS as the results of Refs. 7 and 8. In terms of weighted densities, we obtain →T ←  ξ + ξA 3 16π app → → → n2 → n 2− n 2 − n2 − n 2− n2 + B φ3 = ξ B − n2 3 3  →2 2ξ − ξB ← Tr → n 32 − ξA n2 Tr ← + A n 2 . (20) 3 Condition (iii) amounts to demanding that the kernel in Eq. (16) vanishes when ni → nj for any two i and j. This restriction is necessary to avoid divergences for one and two cavities as well as a quasi-one-dimensional system.7 The same condition is also required to have a stable hard sphere crystal.8 The kernel (17) that is exact in the zero-dimensional limit satisfies a stronger condition than (iii), that is not obeyed by the approximate kernel (19) with ξ B = 0, namely that the kernel vanishes whenever n1 · (n2 × n3 ) = 0. Apparently, violating this condition does not lead to a diminished performance of the functional for hard spheres; however, we will show that it does have an effect for highly anisotropic shapes in this work. Now, we will briefly explore an example system, that of hard parallel cylinders, to determine for which normal vectors this stronger condition is violated by the approximate kernel. In Fig. 1, we depict configurations of normal vectors that we might expect to have a finite contribution to Eq. (16). Configuration (II) has n1 · (n2 × n3 ) = 0; therefore, it contributes to φ 3 without violating the stronger condition. However, the kernel from Eq. (19) is also nonzero for the coplanar triplets of normals (I) and (IIIa) for which n1 = n2 = n3 , which present a violation of the stronger condition, as n1 · (n2 × n3 ) = 0, but not of condition (iii). For all other triplets of normal vectors depicted in Fig. 1, both versions of the kernel vanish irrespective of the value of ξ B . In Sec. III, we will determine whether contributions (I) and (IIIa) present a problem

FIG. 1. Combinations of normal vectors on the surface of parallel cylinders. The illustration projects the outward normal vectors ni of three different bodies into one picture for five representative configurations. The excess free energy ex only shows a correct scaling behavior in case (II), where intersections of two cylindrical shells and one capping disk contribute. In all other cases, the triple scalar product n1 · (n2 × n3 ) accordingly vanishes in the kernel K3 , Eq. (17), of the general third term φ 3 , derived from diapp mensional crossover. The second term in the approximated kernel K3 in Eq. (19), however, depends differently those normal vectors and thus does not vanish in the cases (I-IIIa) without strictly parallel normals.

for more experimentally relevant phases of elongated particles by a scaling argument. III. SCALING ARGUMENTS AND NEMATIC PHASES

In this section we apply the different versions of the functional to verify the scaling behavior for extremely elongated rods at nonzero packing fraction. We compare the emerging constraints to the functional with the results of Sec. II B to obtain an unambiguous expression for the third term of the functional in Sec. III A. Findings for the second term are discussed in Sec. III B in the context of an exact treatment of the ζ correction for HSC. Taking into account all tensorial weighted densities, we obtain the exact mutual excluded volume. The EOS for the parallel system is compared to the result for HC in Sec. III C. Newly gained insights into the functional are applied in Sec. III D to the isotropic–nematic coexistence and interface, originally considered in Refs. 26 and 28 within edFMT. A. Highly aligned spherocylinders

First consider the nematic fluid of HSC with cylinder length L, spherical diameter D, and volume v = π6 D 3 + π4 D 2 L. The corresponding weighted densities can be found in Refs. 27 and 28. The isotropic–nematic transition in the Onsager limit22 of infinite aspect ratio l = L/D occurs at zero packing fraction η = 0. The second virial coefficient, Eq. (7), is sufficient to characterize this transition. The dimensionless concentration c∝ρL2 D remains finite and the nematic phase is only moderately aligned. A more detailed review can be found in Ref. 57. At higher packing η > 0, the scaling of the functional towards l → ∞ needs to be carried out differently to the case η = 0, as then ρLD 2 = ηLD 2 /v remains finite. As a consequence, only orientationally ordered phases are present and higher virial coefficients become important, such that the third term φ 3 in Eq. (5) is not expected to vanish. The goal is to verify, that the free energy  per particle remains finite in these systems, i.e., /ρ = O l 0 . To study the connection between different constituents of the excess free energy from Eq. (5) using at most rank two

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-5

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

(ζ )

tensors, we consider φ2 from Eq. (11) and the most general form  →3  3  3 (ξ ) → ξ1 n2 + ξ2 n2 Tr[← φ3 = (21) n 22 ] + ξ3 Tr ← n 2 16π of the third term, including all possible combinations of scalar → and tensorial weighted densities and neglecting the vector − n2 for the moment. The orientational distribution, which minimizes the density functional from Eq. (1), has the form g(α, cos θ ) =

α exp(−α 2 (1 − cos2 θ )) D(α)

(22)

for the nematic phase, where D(α) denotes Dawson’s integral.26, 27 The excess free energy only depends on the nematic order parameter S(α) = 3 cos2 (θ ) − 1 /2 =

3 3 1 − − 4αD(α) 4α 2 2

(23)

through the tensorial weighted densities.27, 28 The brackets denote an average with respect to the orientational distribution, Eq. (22). The intrinsic order parameter α is defined by the self-consistency equation α2 = −

3 ∂ex (2) 2 = A2 ζ S + A(1) 3 S + A3 S , 2ρ ∂S

(24)

which generalizes  26. The constants   the version from Ref. 1 A2 (η, l) = O l 1 and A(i) with i ∈ {1, 2} arise (η, l) = O l 3 from φ 2 and φ 3 in Eq. (5), respectively, and depend on the explicit form of the functional. Since the nematic order parameter S ≤ 1 from Eq. (23) is bounded, we can expand Eq. (24) in powers of l. We obtain   α 2 =: α12 l + α02 + O l −1 (25) to calculate the expansion of S =1−

  3 −1 6α02 − 3 −2 l + l + O l −3 . 2 4 2α1 4α1

(26)

Using the last two equations, we can efficiently solve Eq. (24) at large aspect ratios l by iteration for the coefficients α k . Now, we insert the scaled nematic order parameter from Eq. (26) into the expressions for φ 3 , given by Eqs. (14) and (15). We find a divergence of the excess free energy   per particle ex /ρ = O l 1 with φ3(T) but not with φ3(TR) , an earlier proposal by Tarazona and Rosenfeld.7 A similar di(ζ ) vergence is observed in φ2 from Eq. (11), when we choose ζ = 2. It is straightforward to show that the ideal gas free energy per particle only has a logarithmic divergence in the orientational entropy ∝ ln α, which does not affect the phase behavior. To generally understand the scaling of the excess term, we extend our analysis to one-dimensionally inhomogeneous systems. As shown in Appendix A, the weighted densities (pa) nν,C (z) of the perfectly aligned cylindrical body, defined in Eq. (A1), constitute the only leading order (LO) terms of the weighted densities. The explicit results can be found in Eq. (B1) of Appendix B, with a general integral C0 (z)O l 1 over an arbitrary one-dimensional density modulation, defined in Eq. (A2). Indeed, we find that ζ = 2 is a necessary

condition or all highly aligned phases, as the LO contribution     1 (ζ ) 2 3 1 −ζ C0 (z)2 + O l −1 (27) φ2 = πρ D 4 8   the second term needs to vanish, such that φ2 /ρ = O l 0 . The third term from Eq. (21) yields     ξ3 3π 2 ρ 3 D 6 ξ2 (ξ ) ξ1 + + C0 (z)3 + O l −1 . (28) φ3 = 16 2 4 The general condition (29) 4ξ1 + 2ξ2 + ξ3 = 0  0 = O l , i.e., a vanishing LO. To preserve the ensures isotropic EOS, which is well represented within the original functional,28 or equivalently the Percus-Yevick EOS in the hard sphere limit, we need to demand (ξ ) φ3 /ρ

9ξ1 + 3ξ2 + ξ3 = 2.

(30)

The necessary requirement to describe the limit l → ∞ or highly aligned systems in general, is given by Eq. (21) and the parameters ξ1 ≡ ξ, ξ2 = 2 − 5ξ, ξ3 = −4 + 6ξ.

(31)

Obviously the original expression of Rosenfeld,2 which has ξ 2 = ξ 3 = 0, does not fulfill this condition. Thus we provide a compelling motivation to use tensorial weighted densities for anisotropic bodies without borrowing expressions from extremely confined systems. However, this argumentation does not provide any information about vectorial contributions. Combination with the general dimensional crossover result from Sec. II B fixes the parameters ξ = 1 as well as ξ B = 0 and ξ A = 3 in Eq. (20). The result ξ B = 0 is consistent with the vanishing LO in Eq. (21). Figure 1 shows that n1 · (n2 × n3 ) = 0 in case (I) of three cylindrical intersections, such that there is no LO contribution in φ 3 when the exact kernel from Eq. (17) is used. The term in the approximated kernel from Eq. (19), which is proportional to ξ B , does not vanish and the free energy becomes divergent for bodies of extreme elongation if we choose ξ B = 0. As a result, we will use Tarazona and Rosenfeld’s7 expression, Eq. (14) instead of the newer version by Tarazona8 which has ξ B = 1 from now on, if not denoted otherwise. Additionally using ζ = 2, we have a nematic excess free energy with the correct scaling dependence. In the limit l → ∞, it exclusively depends on the next-to-leading order (NLO) term  1 3 (32) ⇒ θ = O l− 2 S2 := 3θ 2 /2 = 2α12 l of the order parameter from Eq. (26). Thus, in contrast to a perfectly aligned system, there is still a relation between the density and orientational degrees of freedom. In Appendix A, we define all LO and NLO contributions to the weighted densities of maximal tensor rank two, and minimize the scaled smectic free energy in the limit l → ∞. We present explicit expressions, obtained with a simple cosine modulation, in Appendix B. Noteworthy, the scaling of the angle θ, found  in Eq. (32), is in contradiction to the scaling θ = O l −1 done

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-6

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

in simulations.32 However, in Sec. III B, we provide an approach to higher order terms, which suggests that the scaling behavior of the ζ corrected functional should be incorrect.

contributions gn (γ ), considered in Ref. 27 up to rank 6. The general contribution of rank 2n tensors reads gn =

Πnk cos2k γ ,

(39)

k=0

B. Exact second virial coefficient

We have seen that the LO of the weighted densities may lead to a divergence at nonzero packing fraction for increasing aspect ratio l, corresponding to stronger orientational alignment. In a system of parallel HSC, these effects are even more pronounced, which demands a correct treatment of the functional. We find an exact series representation for the second virial coefficient from tensorial weighted densities in Appendix C. For the nematic phase, there are no hemispherical terms in the combined contribution  [2n] [2n] L2 D = ρ2  , n1 n2 16π n

(33)

of rank 2n tensors. The result of the trivial integrals along the direction of the cylindrical axis is L and the prefactor   1 3  n − 2 2 (34) n = 2π (n + 1) originates from adding the integrals over the polar angle for all tensorial components. The generalized ζ correction from Eq. (12) gives rise to the alternative expression  2 n  (2n − 1)(4n + 1) n − 12 (ζ ) [n] n = m ζ m = π . (35) (n + 1)(n + 1)2 m=1 (ζ )

As we have found closed formulas for n and n , we can explicitly investigate the towards the exact ex∞convergence (ζ ) 2  =  = 4π for summation over pression ∞ n n=1 n  n=1 nt n /(4π 2 ) > 0.9 demands nt ≥ 32 n. An accuracy of n=1 without rescaling and nt ≥ 6 with the generalized ζ correction, when the expansion is truncated at n = nt . We further find from the comparison to 1 = 2π 2 , that the corrected rank two functional with ζ = 2 is exact for the nematic phase of parallel HSC. Note, that this is not necessarily the case for more complex phases with additional hemispherical contributions, see, e.g., the perfectly aligned smectic phase in Sec. IV B. The complete nematic EOS βpv = y + 3y 2 + 3y 3 , y =

n 

η 1−η

(37)

of two HSC. We find the difference vexcl (0) − vFMT (0) π = (2 − ζ ) 2L2 D 8

which generalizes Eq. (24). For the nematic order parameter S → 1, some divergence of the derivative of ζ (S) seems to be found when it is used to minimize the error in the excluded volume. However, the result is only known numerically, while the parabolic approximation ζpar (S) = (5 + 3S)/4 has a finite derivative at S = 1 and thus behaves similar to ζ = 2.27 Given the more general form ζ  (S) ∝ (1 − S)−1/2 as S → 1, such that (1 − S)−1 ∝ S2−1 ∝ θ 2 −1 ∝ α 2 √ ∝ A2 ζ  (S) ∝ l/ 1 − S as l → ∞, S → 1 would imply the correct scaling θ 2 ∝ l−2 , but still contain the wrong order parameter. To improve the analytic approximation for the numerical data of ζ (S)27 beyond ζpar (S), we therefore suggest a √ more general expression which includes a term ∝ 1 − S. C. Scaling of parallel hard cylinders and cubes

(36)

differs only in the last term from a result for parallel HC.15 The finding of ζ = 2 to be exact for parallel HSC is in perfect agreement with an argumentation related to the excluded volume 4 vexcl (γ ) = 2L2 D sin γ + 2π LD 2 + π D 3 3

with the constants Πnk specified in Appendix D. We see from summation over n ,that the extended deconvolution of Eq. (10) is indeed capable to describe the exact second virial coefficient. Re-examining the scaling of the orientational degrees of freedom in Sec. III A, we expect from comparison of the orientation dependence in Eqs. (37) and (39), that an exact representation of φ 2 as in Eq.  would give rise to the  (8) correct scaling behavior θ = O l −1 . This can be seen from Eq. (32), as a different orientation-dependent NLO contribution proportional to θ needs to be considered in the case of an exact second term. However, in case of such an exactly included second virial coefficient, a similar calculation as shown in Appendix A, reveals that the third term would lose its orientation dependence via θ 2 in the limit l → ∞. Interestingly, in Ref. 27 a semi-empirical running correction ζ (S) of φ2[2] was proposed which, in a different way, also could restore the correct scaling of θ . This is due to the selfconsistency condition   1  (2) 2 2 2 (40) α = A2 ζ (S)S + ζ (S)S + A(1) 3 S + A3 S , 2

(38)

between the exact result the tensorial contribution from Eq. (53) in Ref. 27 with the intermolecular angle γ = 0. This result can be reproduced with the full series over all tensorial

The hard cylinder (HC) fluid is related to HSC by a common limit for l → ∞, as illustrated in Fig. 2. A HSC fluid yields the well-studied hard sphere limit for l → 0. HC have a richer phase behavior, as their shape interpolates between extremely oblate disks and prolate rods in the respective limits l → 0 and l → ∞. The most isotropic shape, obtained for l = 1, is still quite complex due to elliptic integrals in the mutual excluded volume of a HC.22 The parallel HC system, however, possesses a universal scaling behavior along the cylindrical axis, which makes the system independent of the aspect ratio l. The phase behavior can be accessed from the perspective of two-dimensional FMT functionals,2, 7 scaled to a three-dimensional version for parallel HC.14, 15 It is logical that the EOS obtained from edFMT should be consistent with this alternative approach, as both go back to the first principles

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-7

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

for the parallel HC system from edFMT. This EOS is to be compared to the result

HSC

sphere

l→0

βpv = y + 3y 2 + 2y 3

l→∞

l=0

l=∞

line disk HC FIG. 2. Scaling behavior of a hard spherocylinder (HSC) with cylinder length L and diameter D. At vanishing aspect ratio l = L/D = 0, a HSC reduces to a sphere. The phase behavior of parallel hard cylinders (HC) is independent of the aspect ratio and the limit l → ∞ of both a HC and a HSC is an infinitely thin rod. Hence, a perfectly aligned HSC can be scaled towards thin platelets, i.e., a HC with l → 0.

of FMT, namely scaled particle theory.4 Likewise, edFMT is constructed to reduce to the Percus-Yevick3 EOS in the hard sphere limit. Thus, a study of HC can reveal additional insight into the functional. Different to a HSC, considered in Sec. III A, the third (ξ ) term φ3 from Eq. (21) has to obey two different scaling laws for extreme anisotropy of a HC, such that the free energy per particle /ρ in the parallel system does not depend on the aspect ratio. Therefore, not only the LO, which is equivalent to the HSC system, but rather all terms in , which are not proportional to the inverse aspect ratio l−1 , have to vanish. Using this additional condition, we can confirm that the remaining parameter in Eq. (31) should be ξ = 1, without considering dimensional crossover in Sec. II B. The drawing in Fig. 1 explains that both divergences for l → ∞ and l → 0 are equivalent to over-weighted contributions of respective configurations (I) and (IIIa) in the approximated kernel from Eq. (19). Now that we have two independent arguments for using third term from Eq. (14) instead of Eq. (15), we can study the parallel HC system consistently. We use the weighted densities of a HC fluid  π π  n3 = ρ LD 2 = η , n2 = ρ π LD + D 2 , 4 2   D L +π , n0 = ρ, n1 = ρ 4 8 π  π → → LD(2 + S) + D 2 (1 − S) , n 2 )22 = ρ (← n 2 )11 = (← 6 6 (41) π  π → LD(1 − S) + D 2 (2S + 1) , (← n 2 )33 = ρ 3 6   π 1 1 ← → ← → ← → LS − DS , ( n 1 )11 = ( n 1 )22 = − ( n 1 )33 = ρ 2 8 32 with S = 1 and φ2(ζ ) from Eq. (11) with ζ = 2 to obtain the nematic EOS 3 9 βpv = y + (2 + π )y 2 + y 3 4 2

(42)

(43)

from Ref. 15. This reveals two severe drawbacks of the present functional. In contrast to the HSC system, which yields the EOS presented in Eq. (36) for l → ∞, the ζ correction does not appear to be consistent, although the value ζ = 2 unambiguously guarantees the correct scaling. This can be understood by the more advanced expression for the excluded volume of HC.22 The different orientation-dependent terms cannot all be restored correctly by one single prefactor. In addition, a comparison of the third terms in in Eqs. (36) and (42), reveals that a HSC at l → ∞ is not identically represented as a HC, which disagrees with intuition and the scaling in simulations.32 Interestingly, we can reproduce Eq. (43) with HSC, choosing ξ = 0 in Eq. (31), although this third term has the wrong dimensional crossover properties and the HC system would not scale properly towards l → 0. Recall that, as long as we only consider rank two tensors, the third term is unambiguously determined, except for the prefactor, which is necessary to retain the hard sphere limit. Specifically for parallel HC, we find a deviation by a factor of 4/9 between the last term in Eq. (42) and the scaled particle result, Eq. (43). A similar discrepancy for parallel hard cubes with rounded edges was remedied by the introduction of a semiempirical prefactor to the third term of the original edFMT functional, interpolating between the third virial coefficients of cubes and spheres.30 Considering parallel hard cubes with our new expression for the third term, Eq. (14), we obtain the same incorrect result as from Eq. (15). The full kernel in Eq. (17) yields the correct contribution to the scaled particle EOS from Ref. 13. Note, that for this system there are no lost cases in the third term. However, the expression for the kernel is rather complicated and thus impractical for orientational disordered systems or continuously curved bodies. To this end, we discuss the empirical modification φ3(TR) → 49 φ3(TR) for hard thin platelets at the end of Sec. III D.

D. Results at isotropic–nematic coexistence

In Sec. III A, we have seen that the variation of the prefactors in Eq. (21) has a significant impact on the validity of the functional for extremely anisotropic bodies. After imposing the condition from Eq. (30), the EOS of homogeneous and isotropic bulk fluids does not depend on the remaining parameters and matches the result of FMT in the fundamental hard sphere limit. For now, we only consider systems with moderate orientational order, where the discovered divergence does not occur. Still, the nematic phase may be strongly affected when the dependence on the order parameter in the third term φ 3 is modified. The isotropic–nematic transition of infinitely long rods at η = 0 and S < 1 is well described by Onsager’s low-density expression and thus independent of the third term. For decreasing aspect ratio, the transition is shifted towards higher densities and this term becomes more important for the phase behavior.

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

0.8

0.5

ηI , ηN

0.4

0.45

ζ=1.6, Eq. (15) ζ=5/4 ζ=2 simulation

0.4

0.35

0.3

η

0.6

L/D=5

0.4

SN

064103-8

nematic

ηI , ηN

isotropic

0.1

0.15 D/L

0.2

0.12

0.25

FIG. 3. Isotropic–nematic coexisting densities ηI and ηN of HSC, based on the edFMT functional with the third term from Eq. (14). We show results as a function of the inverse aspect ratio l−1 = D/L for ζ = 5/4, which yields the smallest quadratic error for the excluded volume, and ζ = 2, which is exact for parallel rods, with dotted and solid lines, respectively. Comparison is made to computer simulations and the original results of Ref. 26 using ζ = 1.6 and Eq. (15). The difference between coexisting densities and the behavior for higher aspect ratios is significantly improved with ζ = 5/4 as emphasized in the inset.

We see in Fig. 3 that the results obtained with the third term φ3(TR) from Eq. (14) do not match the simulations32 as well as the ζ = 1.6 corrected original version.26 However, there is no justification for the latter empirical choice other than this nice fit to simulations, while there is a wellmotivated value ζ = 5/4, which minimizes the quadratic difference in the mutual excluded volume.26, 27 As already mentioned in Ref. 28, there is another motivation to choose ζ = 5/4 as the first order result of a systematic expansion of the Onsager functional,58 i.e., the generalized ζ correction defined in Eq. (13). We use this expansion to show in Fig. 4, that the ζ = 5/4 correction is indeed the best approximation to describe the isotropic–nematic transition including at most rank two tensors. The addition of higher order terms up to maximal tensor rank nt systematically increases the difference δηIN = ηN − ηI of coexisting densities and the coexisting nematic order parameter SN . As also shown in Ref. 58, we may consider the presented results at isotropic–nematic transition of HSC with nt = 10 as converged, although we need to include more terms for higher ordered problems. We find the most accurate description when we use the third term φ3(TR) . On the other hand, it needs to be pointed out that the higher order results clearly deviate from the empirical fit using ζ = 1.6, when we use the original version of edFMT, which includes φ3(T) from Eq. (14). We generally see that larger values of ζ for rank two tensors induce lower transition densities. Most notably, the value ζ = 2, which is excellent for highly aligned systems shown in Sec. III A, yields a poor description of the isotropic–nematic transition in Fig. 3. Note, that the running ζ (S) correction proposed in Ref. 27 interpolates between the two presented semi-empirical choices ζ = 5/4 and ζ = 2. Now that we modified the functional, this approximation turns out to be more appropriate than in combination with φ3(T) . Figure 4 reveals, that the results with the numeric values for ζ (S) are relatively close to the higher order results for both choices of φ 3 . The parabolic approximation ζpar (S), on the other hand, produces

ζ

2

4 6 2nt

ζ

8 10

2

4 6 2nt

8 10

0.4

FIG. 4. Isotropic–nematic coexisting densities ηI and ηN (left) and coexisting nematic order parameter SN (right) for the aspect ratios l = 5 and ζ [n ]

l = 20 of HSC. The second term φ2 t from Eq. (12) is truncated after different orders r = 2nt of the generalized ζ correction, defined in Eq. (13) and (TR) applied together with φ3 from Eq. (14) (circles) and Tarazona’s original (T) term φ3 from Eq. (15) (crosses). We also show second order results with (ζ ) (TR) (squares) φ2 from Eq. (11) using the numeric ζ (S) correction27 for φ3 (T) and φ3 (plus sign). Diamonds denote the original results for the latter expression using the constant parameter ζ = 1.6.26 The generalized correction of rank two tensors is equivalent to ζ = 5/4. We may deduce from the figure that an expansion up to tensors of rank 10 is sufficient for an adequate description of the isotropic–nematic transition of HSC. We compare the DFT results to simulation data32 and Parsons-Lee theory43, 44 shown with the thick and thin dotted lines, respectively.

results similar to ζ = 1.6 and seems to overestimate the width of the coexistence.27 Note that the implementation of a local, orientation-dependent correction, used to correct the absence of non-local weighted densities, would be inconsistent for problems with varying order parameter.27 Within the restriction to rank two tensors, we already see in Fig. 5, that φ3(TR) improves the description of the isotropic–nematic transition on the level of the coexisting density difference δηIN . We also show the interfacial tension γIN [ρ] = ([ρ] + pV )/A, minimized by the width δ of the interface and the offset z of the order parameter profile. The

0.15

ζ=1.6, Eq. (15) ζ=5/4

simulation DFT Onsager

0.8 0.6

0.1

l δηIN

0.05

0.6

0.14

L/D=20

βγIN (L+D)D

0.1 0 0

0.8

0.16

SN

0.2

0.4 0.05

0.2 0

0

0.1

D/L

0.2

0

0.1

0.2

D/L

FIG. 5. Isotropic–nematic interfacial tension γIN for different aspect ratios of hard spherocylinders and related differences δηIN = ηN − ηI of coexisting densities, which are shown in Fig. 3. Our new results, using Eq. (14) and ζ = 5/4 (circles), show strongly increased values for γIN , compared to the original values from Ref. 28 (diamonds) at small aspect ratios. As remarked in Ref. 28, the value of the ζ parameter is of minor importance, which becomes manifest in the similar Onsager limit. Comparison is made to grand-canonical Monte Carlo simulations59–61 (downward triangles), the DFT by Somoza and Tarazona42, 62 (upward triangles), and to Onsager theory63 (squares).

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-9

Wittmann, Marechal, and Mecke

associated density modulation28   ρ (z,  ) = ρN h(z) + ρI (1 − h(z)) g (α(z), sin θ cos φ) (44) √ with α(z) = αN h(z − z), h(z) = (1 + tanh (z/δ))/2 and the orientational distribution from Eq. (22). The equilibrium nematic director is likewise parallel to the interface, described by the substituted argument in g(α, cos θ ), which yields the maximum at θ = π /2. The comparison to δηIN for different aspect ratios confirms the scaling relation derived in Ref. 28 for the interfacial tension but the result βγIN (L + D)D = 0.0654 at l = 15 is still smaller than βγIN (L + D)D = 0.1 from simulations.59–61 The interface width δ = 0.80L compares well with the simulation result δ = 0.71L,61 improving upon the original approach.28 The comparison with the results of the Somoza-Tarazona DFT62 also indicates clear improvement with the current edFMT at low aspect ratios L/D. The deviation to the result in the Onsager limit,63 accounts for the failure of the ζ correction to approximate the second virial coefficient, as the functional does not depend on φ 3 in this limit. Note, that this correction is exclusively based on the mutual excluded volume. Therefore, it could be expected to be even more reliable for highly elongated HSC than for short rods, in particular with ζ = 5/4. In fact, the value of the ζ parameter barely has any influence on the value of the interfacial tension for all aspect ratios.28 We further find the same qualitatively correct behavior of the interfacial tension as a function of the director orientation with φ3(TR) , when we use a more general orientational distribution in Eq. (44).28 Finally, we turn to the opposite end of the shape spectrum and consider a fluid of infinitely thin platelets, i.e., HC of vanishing aspect ratio. In contrast to long rods, the isotropic– nematic transition of disks depends sensitively on higher virial coefficients. As a result, three body overlaps, as generally described by φ 3 in Eq. (16), are important in this case. This system is also of particular interest from the viewpoint of FMT, as there exists a fundamental measure functional specifically derived for such a particle shape.21 The equivalent expressions to φ 2 and φ 3 from Eq. (5) in that functional contain mixed weight functions, which depend on more than one body. Within edFMT, we use the HC weighted densities in Eq. (41) with L = 0. The excluded volume vexcl = π/2D 3 sin γ of platelets22 validates the generalized ζ correction to study this system, which is ideal to test empirical modifications of Eq. (14). Indeed, Table I reveals that the third term (4/9)φ3(TR) , adapted to match the EOS of parallel HC,15 also yields the best description of the thin platelet isotropic– nematic transition, when compared to simulations64, 65 and DFT.21 Interestingly, both Rosenfeld-like versions φ3(T) and φ3(TR) , respectively, overestimate and underestimate the difference of coexisting densities and the nematic order parameter significantly. Using the generalized ζ correction within Onsager’s second virial approximation, which only considers Eq. (7) for the functional, we find slightly different results than those stated in Refs. 21 and 64. In future work, we need to find a sensible modification of the third to improve the description of more general systems. Unfortunately, it ap-

J. Chem. Phys. 141, 064103 (2014) TABLE I. Isotropic–nematic coexisting densities σ = ρD3 and nematic order parameter SN for hard thin platelets with different versions of edFMT, compared to simulations64, 65 and DFT results.21 As for spherocylinders, the second term is corrected by ζ = 1.6 for the original third term from Eq. (15) (TR) and ζ = 5/4 with Eq. (14). Additionally, the correction factor 4/9 to φ3 is chosen, such that the EOS from Ref. 15 is recovered for parallel platelets (cylinders). More accurate results, also for the Onsager approximation, are calculated with the generalized ζ correction up to tensors of rank 16. Work Onsager, nt = 8 (T) φ3 , ζ = 1.6 (T)

φ3 , n t = 8 (TR) φ3 , ζ = 5/4 (TR) φ3 , nt = 8 4 (TR) , ζ = 5/4 9 φ3 4 (TR) , nt = 8 9 φ3 Frenkel64

Marechal65 Esztermann21

σI

σN

SN

5.334 3.451

6.791 4.496

0.792 0.803

3.737

5.508

0.900 0.171

2.948

2.991

2.944

3.009

0.211

3.771

3.917

0.302

3.747 3.78 4.03 3.344

4.017 4.07 4.34 3.680

0.407 0.37 0.52 0.492

pears very unlikely, that a general, less empirical approach can maintain computational efficient structure of FMT. IV. SMECTIC PHASES OF SPHEROCYLINDERS

Now that we can account for the correct scaling behavior, we present consistent results for smectic phases, using a generalized FMT functional, for the first time. Section IV A introduces the smectic modulation function and calculational methods. The analytic analysis of perfectly aligned HSC in Sec. IV B emphasizes the unambiguous need for our newly proposed version of edFMT. Finally, in Sec. IV C, we present the results for the phase diagram and EOS for HSC including the smectic-A phase with the appropriate functionals. A. Approaches to smectics

Characteristic for smectic phases, is an inhomogeneous spatial density profile with periodicity d = dL l and amplitude λ. The average over one smectic layer gives rise to the average free energy density  1 d = dz(z). (45) d 0 The nematic phase is obtained for λ = 0. For our numerical calculations, we will use the spatial density profile42   2π z (46) ρ(z) ∝ exp λ cos d d in Eq. (2), normalized such that ρ = 0 dzρ(z)/d. The orientation dependence is taken as g(α, cos θ ) from Eq. (22), which minimizes the nematic free energy. The function has its maximum for θ = 0 and describes a director perpendicular to the smectic layers, which is characteristic for smecticA alignment. The intrinsic order parameter α may additionally depend on z by a similar modulation as the density in Eq. (46). However, we will stick to the orientational decoupling approximation and assume that α is a free parameter,

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-10

Wittmann, Marechal, and Mecke

which does not depend on z. The scaled orientational distribution g2 (α, θ ) in the limit l → ∞ can be directly derived from a variation of the smectic free energy, as carried out in Appendix A. With the simplistic density profile 2π z =: ρ + ρ (z) (47) d for the smectic modulation, which is the expansion of Eq. (46) up to first order in λ, the weighted densities nν of highly aligned HSC can be obtained analytically. The explicit expressions can be found in Appendix B. The resulting free energy density (z) may be expanded for small amplitudes λ to analytically calculate the average     1 d k dz 0 + k (z)λ (48) = d 0 k=2 ρ(z) = ρ + ρλ cos

over one smectic period afterwards. The zero order 0 is equal to the nematic free energy. Hence, we can perform a stability analysis of the smectic free energy given by Eq. (48). A more elaborate analysis, as performed in Refs. 34, 35, 37, and 52, involves additional Fourier components in Eq. (47). Alternatively, we can determine the order of the transition from an explicit numerical evaluation of the weighted densities using Eq. (46). To do so, we minimize the smectic free energy at given packing fraction η with respect to the parameters λ, d, and α of the trial profile from Eq. (46) with the orientational distribution in Eq. (22). The smectic grand potential corresponds to the minimum of Eq. (1) at given chemical potential μ. To determine the phase coexistence, we simply need to find the value, which is equal to the nematic grand potential at the same μ. Before presenting the results of this computational intensive calculation, we focus on parallel HSC, which are a good analytic reference system to verify the behavior of the orientational disordered system for specific approaches. B. Perfectly aligned system

In a system of perfectly aligned rods, there is only one homogeneous phase, which draws special attention on the nematic–smectic transition. Using only the weighted densities of parallel HSC from Appendix B, obtained with the profile from Eq. (47), the full smectic free energy can be calculated analytically for arbitrary aspect ratio. Thus, it is convenient to study smectics in the perfectly aligned system first and quantify the influence of the ζ parameter and different choices for φ3. The simplified system further gives rise to a general calculation of higher orders in the tensor expansion, Eq. (8), for the smectic phase. In Appendix C, we find the contributions  [2n] [2n] D3  C (z)(C0 (z) + Hn (z)) (49) n1 n2 (z) = ρ 2 16π n 0 of rank 2n tensors to the exact second term in Eq. (8). Explicit expressions for the contributions of cylindrical C0 and hemispherical Hn integrals, involving the density modulation from Eq. (47), are given in Appendix B. However, the contracted tensorial contributions in Eq. (49) can only be evaluated when we neglect Hn , which vanishes for the nematic phase, studied equivalently in Sec. III B. As shown in Appendix B, the series

J. Chem. Phys. 141, 064103 (2014)

involving Hn appears to converge for the smectic phase, but cannot be evaluated. For a mathematically convenient approach to the exact hemispherical contribution H in Appendix B, the sum over n needs to be evaluated before the spatial integration, i.e., before are calculated. In this sense we the weighted densities n[2n] ν would have to define a mixed weighted density  ∞   (1)  (2)  ω[r] (r ) ω[r] (r ) . dr dr ρ(r − r)ρ(r − r) n12 = r=2

(50) In Appendix C, we show, that we can write π D3 C (C + H) (51) 4 0 0 without facing difficulties arising from an alternating series as in Eq. (49). Recall that we find 1 = 2π 2 from Eq. (34), which means that the smectic phase with 2H1 = H is not exactly described by using ζ = 2 for n = 1 in Eq. (49), unlike the nematic phase with H1 |λ=0 = H|λ=0 = 0. We show in Appendix B that both hemispherical terms scale out in the limit l → ∞, such that the parameter ζ = 2 yields the exact smectic free energy. For all considered variants of edFMT, the Landau expansion from Eq. (48) is invariant under  transition λ → −λ  the and becomes  = 0 + 2 λ2 + O λ4 . We solve the bifurcation condition 2 = 0 with the free parameter d, to determine the packing fraction η of the second order transition, as 4 > 0 can be verified by an explicit calculation. As argued in Appendix A, orientational disordered HSC of infinite length are almost perfectly aligned and can be treated similarly. The orientation-dependent contributions to the weighted densities, which become relevant in this case, are also given in Appendix B. The resulting expansion of the free energy only contains even terms and we find a second order nematic– smectic-A transition, even when we use higher orders of the Taylor expanded Eq. (46) for the density modulation. The nematic–smectic transition of parallel HSC, shown in Fig. 6, confirms both the qualitative simulation results51 of a second order transition with decreasing density for higher aspect ratio, and the scaling behavior for highly aligned systems considered in Sec. III A. The divergence of the transition density for l → ∞ due to the original third term φ3(T) (ζ ) from Eq. (15) or a correction of φ2 from Eq. (11) with ζ = 2 becomes manifest at a quite small aspect ratio. The only qualitatively correct truncation result, using ζ = 2, is in remarkable agreement with the evaluation of the infinite series, φ2[∞] from Eq. (8). However, the inset of Fig. 6 shows a very small non-monotonic behavior due to the approximated term. At l → ∞, we find η = 0.3024 in perfect agreement between those two approaches, as expected in Sec. IV B. In general, the transition density is significantly smaller than observed in the simulation,51 as shown quantitatively in Table II. One reason is that, for uniform phases, FMT reproduces the result of scaled particle theory, which is known to overestimate the fluid pressure. Therefore, the functional for parallel HC15 locates the nematic–smectic transition near edFMT for long HSC, although the functionals differ in their third term, which is manifest, e.g., in Eqs. (36) and (43). The White Bear mark II n12 = ρ 2

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-11

Wittmann, Marechal, and Mecke

0.6

J. Chem. Phys. 141, 064103 (2014)

3

ζ=2, Eq. (15) ζ=5/4 ζ=2 n12 simulation

smectic

Eq. (15)

2

η

η(z)

0.5

ζ=2

2.5

ζ = 5/4

1.5

ζ=2

1

0.4

0.5 0.3 nematic

0.2 0

0.5

0.3024 0 1 D/L

0 0 0.1 1.5

version of edFMT,56 which reduces to the more accurate Carnahan-Starling EOS54 in the hard sphere limit, and even the general third term given by Eqs. (21) and (31), with ξ as a free parameter (spoiling the dimensional crossover), do not change the location of this transition significantly. Remarkably, the edFMT predictions, 2.18 > dL > 1.30 for 1 < l < ∞, for the smectic layer spacing d = dL L deviate by less than 0.1 from simulation results.49, 50 In order to further motivate the choice ζ = 2 and φ3 = φ3(TR) , we consider briefly a system of parallel HSC confined between two hard walls with separation H. Due to the hard walls at z = 0 and z = H, the z-component zi of the center-ofmass of a spherocylinder i satisfies (L + D)/2 < zi < H − (L + D)/2, so the Mayer f-function for two HSC with coordinates (xi , yi , zi ) is equal to the Mayer function of two hard disks with coordinates (xi , yi ). In the z-direction, the particles behave like an ideal gas, so that the density profile η(z) := vρ(z) should TABLE II. Location ηNS and layer spacing d = dL L of the nematic–smectic transition of parallel spherocylinders at l = 5 and in the limit l → ∞ (last (TR) two columns). The edFMT values are obtained with the third term φ3 from Eq. (14). For the second term, we use the two choices ζ = 5/4 and ζ = 2 to correct missing terms in the expansion from Eq. (8). The motivation for the latter value is a non-divergent smectic free energy. We also include the exact second term with n12 from Eq. (50) and the White Bear mark II structure of the functional.56 The results are compared to simulations,50, 51 DFT,36, 39, 40 and the fundamental measure functional for hard cylinders.15 (∞)

(∞)

ηNS

dL

ηNS

dL

ζ = 5/4 ζ =2 n12 WBII, ζ = 2

0.373 0.304 0.305 0.305

1.374 1.475 1.472 1.478

∞ 0.302 0.302 0.306

... 1.299 1.299 1.302

Stroobants50 Veerman51 Somoza39, 40 Poniewierski36 Capitán15

0.41 0.41 0.35a 0.28a 0.314

1.39a 1.39a 1.48a 1.59a 1.302

0.35 0.44 0.36a 0.28 0.314

1.27

a

Valued are read off from graphs.

2

3

4

5 z/D

6

7

8

9

10

2

FIG. 6. Nematic–smectic transition density for parallel hard spherocylinders. The transition line diverges for l → ∞, when we use Eq. (15) for the third term φ 3 and ζ = 2. This is shown exemplarily with dotted-dashed and dotted lines, respectively. The dashed lines are calculated with Eq. (14), corrected by ζ = 2, which becomes exact at l → ∞ but shows a non-monotonic behavior revealed in the inset. The exact second term, evaluated via Eq. (50), yields a very similar but monotonous behavior. However, significant deviations to simulation results51 can be seen.

Work

1

1.30a 1.39 1.302

FIG. 7. The dimensionless density profile η(z) for parallel hard cylinders as a function of the distance z from the bottom wall, obtained from different versions of the functional. The average packing fraction equals η = 0.3, the aspect ratio is l = 5, and the wall separation equals H = 10D. For these parameters, the exact solution (see text) for the density profile is η(z) = 0.75 for 3 < z/D < 7 and η(z) = 0 otherwise. The lines belong to the same functionals as in Fig. 6.

be that of a single broad layer with  ηH L+D 40, while we do not find any stable smectic phase at a sensible density at all for l > 2 with the original third term φ3(T) . Note, that a positive slope for the transition density as a function of the aspect ratio is also observed in Refs. 35, 37, 38, and 47 but turns out to be wrong as clarified in simulations.32 For both values of the ζ parameter used here, we see a change of the phase transition from first to second order with increasing aspect ratio. This behavior had been assumed to be true for a long time, although predictions of the tricritical point reach from l = 5.938 to l = 50.41 Despite this conjecture, Bolhuis and Frenkel found hints for a first order transition at l → ∞, which was finally observed in a succeeding simulation.33 Unfortunately, we clearly obtain a second order transition in this limit, using the exactly scaled orientationdependent weighted densities in Appendix B. The first order transition, which we generally observe at small aspect ratios, shows a better quantitative agreement with the simulation results32 and also the Somoza-Tarazona functional42 when we choose ζ = 5/4. Note that there is no reason to assume that the approximation with ζ = 2 yields results that are as accurate as for the perfectly aligned system studied in Sec. IV B. The EOS for aspect ratio l = 5, shown in Fig. 9, gives a similar insight. The curves for ζ = 5/4 compare nicely to simulations,31 although the regions of the phase transition are different. As predicted in Sec. III D, this becomes worse for ζ = 2. In addition, we see that the nematic pressure from the original expression for φ 3 with the empirical ζ = 1.6 is clearly too high, although the location of the isotropic– nematic transition fits perfectly.26 Hence, the effect of different choices for φ 3 in Eq. (5) (and the ζ parameter) clearly dominates the shortcomings of scaled particle theory, as the isotropic pressure is only slightly too high,28 compared to the deviations in the nematic branch. The characteristics of the smectic branch, which is only stable with φ3TR , are independent of the ζ correction. In contrast to the isotropic– nematic transition, we find that the smectic-A phase becomes

FIG. 9. Equation of state of hard spherocylinders with aspect ratio l = 5. The edFMT results are denoted by lines with the same labels as in Fig. 3. In addition, we include the isotropic and nematic branch from the original func(T) tional, using φ3 from Eq. (15) and ζ = 1.6. Horizontal lines denote coexisting pressures and separate stable from metastable branches. The simulation results from Ref. 31 are displayed with different symbols, which correspond to individual phases.

metastable above the nematic bifurcation point. Increasing the density, we find a larger amplitude λ (not shown), which accounts for stronger layering effects, whereas the layer spacing d shows a non-monotonic behavior for l = 5. Table III reveals that, at isotropic–smectic-A and later nematic–smectic-A coexistence, both parameters decrease with the aspect ratio, as it was observed with the Somoza-Tarazona functional.41 In Table IV, we compare the coexisting densities and smectic layer spacing for l = 5 and l → ∞ to simulations32, 33 and different DFT results.34, 38, 41, 42 Showing the HSC phase diagram in Fig. 10, we can conclude on the limitations of the current approximations in the (ζ ) second term φ2 . The only value ζ = 2, which is qualitatively correct for infinitely long rods, suffers from quantitative shortcomings at finite aspect ratio. The location of the isotropic– nematic–smectic-A triple point at l = 2.90 and l = 2.95 for ζ = 5/4 appears to be quite universal and is similar to the simulation prediction32 and other DFT results.38, 41, 45 Apparently, from Figs. 8 and 10, the correction with ζ = 5/4 starts to fail for the smectic-A phase above the tricritical point at 7 < l TABLE III. Densities and order parameters at the isotropic– (for SI = 0) and nematic–smectic-A transition for different aspect ratios l of hard spherocylinders. We use the edFMT approaches described Table II. We also present smectic periodicity d = dL L and amplitude λ in the smectic modulation profile. For l → ∞ the transition density and parameters are determined analytically, using scaled weighted densities from Appendix B. l

ζ

ηI/N

SI/N

ηS

SS

dL

λ

2 3 5 10 20

5/4 5/4 5/4 5/4 5/4

0.487 0.464 0.459 0.459 0.485

0 0.642 0.869 0.946 0.979

0.564 0.522 0.481 0.459 0.485

0.915 0.914 0.921 0.946 0.979

1.47 1.34 1.25 1.17 1.12

10.5 7.16 2.27 0.00 0.00

2 3 5 10 20 ∞

2 2 2 2 2 2

0.444 0.415 0.400 0.382 0.370 0.360

0 0.644 0.858 0.935 0.968 1

0.537 0.485 0.435 0.389 0.370 0.360

0.916 0.916 0.924 0.942 0.968 1

1.50 1.38 1.30 1.25 1.22 1.20

8.66 5.32 2.42 0.89 0.00 0.00

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-13

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

TABLE IV. Coexisting densities ηN/S order parameters SN/S and smectic layer spacing d = dL L for hard spherocylinders of aspect ratio l = 5 (∞) and transition density ηNS for infinitely long rods, obtained with the same edFMT version as in Table II. We compare the results to DFT34, 38, 41, 42 and simulations.32, 33 Note, that Poniewierski used the virial expansion in the limit l → ∞,34 as the density functional from Ref. 38 predicts an increasing transition density. Work

(∞)

(∞)

ηN

SN

ηS

SS

dL

ηNS

dL

ζ = 5/4 ζ =2

0.459 0.400

0.869 0.858

0.481 0.435

0.921 0.924

1.25 1.30

∞ 0.360

... 1.20

Frenkel32, 33 Poniewierski34, 38 Somoza41 Velasco42

0.453 0.526 0.507 0.477

0.74 0.875 0.912

0.482 0.549 0.554 0.503

0.90 0.898 0.985

1.41a 1.26a

0.426 0.417 0.393

1.24 1.21a

a

Values are read off from graphs.

< 8. It is noteworthy that, for ζ = 2 the tricritical point 12 < l < 13 is at higher aspect ratio, which also accounts for a better approximation of the second virial coefficient for increasing alignment. Note, that the orientational decoupling approximation becomes more reliable for higher aspect ratios,53 such that our prediction of a tricritical point is not a result of this assumption. V. DISCUSSION

In the present paper, we demonstrated that the expression for the third term given by Eq. (14) is the only combination of rank two tensors which can provide the correct scaling of the functional for extremely anisotropic hard bodies and simultaneously account for the right dimensional crossover. The implementation leads to a significantly improved phase behavior of HSC, compared to the previously used expression, Eq. (15), which is more consistent for hard spheres. In particular, the isotropic–nematic transition is now well described without a purely empirical correction and the interfacial tension nearly doubles its value for small aspect ratios approach-

smectic-A

0.5

η

0.4 0.3

isotropic

nematic

ζ=5/4 ζ=2 simulation I-N simulation N-Sm

0.2 0.1 3

5

10 l =L/D

20

40

FIG. 10. Comparison of the hard spherocylinder phase diagram, obtained with two distinct values of the ζ parameter for the modified edFMT functional with Eq. (14), in comparison to simulation results.32 Lines and symbols are equivalent to respective results for the isotropic–nematic and nematic– smectic-A transitions from Figs. 3 and 8. For small aspect ratios l < 2.95 and l < 2.90 we find a first order isotropic–smectic-A transition, respectively with ζ = 5/4 and ζ = 2. The conflict between these approximations and possible improvements are discussed in the text.

ing the results of other works. We also show, that the new expression enables the study of smectic phases, which are unstable otherwise. In this context, a conflict between two well-motivated choices for the ζ parameter has emerged: On the one hand, the established, semi-empirical value ζ = 5/4 provides good results for dilute nematic phases and short rods, which are not strongly aligned. On the other hand, we found from the same scaling argument as for the third term that we need to require ζ = 2 to accurately describe the increasing orientational alignment at high aspect ratio. This correction is exact for the artificial systems of perfectly aligned nematics and smectics of infinitely long rods, which can be shown from an analytic evaluation of higher order terms. However, for differently shaped bodies, even as simple as a HC, such a correction is not sufficient to provide the exact result. The interfacial tension for l → ∞ and qualitative failures in the description of the nematic–smectic-A transition, such as a tricritical point at moderate aspect ratios, exhibit the drawbacks of the ζ correction even for HSC. We need to improve the representation of the second virial coefficient to draw a completely consistent phase diagram. The running ζ (S) correction could interpolate between a quantitative description of short rods and the convergence towards a finite transition density at increasing aspect ratio. Therefore, it is doubtful that the location of the tricritical point will change significantly. Higher rank tensors in the second term instead should subsequently shift the tricritical point to higher aspect ratios. The expected slow convergence can be improved with the generalized ζ correction from Eq. (12) as the neglected contribution of odd tensors should be small for symmetric bodies. A more consistent alternative, the expansion of Iijκ from Eq. (10) into spherical harmonics,67 includes these odd terms, and may be readily applied to other particle shapes. However, the convergence of both expansions, which suggest ζ = 5/4 at first order, is still relatively slow for high order parameters, which leads to a large computational effort. Furthermore, only purely empirical corrections, such as ζ = 2 for rank two tensors, could completely overcome the divergence at l → ∞. Therefore, we postpone the implementation of additional terms. It is desirable to derive an exact representation of the second virial coefficient, as in Eq. (50) for parallel HSC. Finding a generalized expression for such a mixed weighted density will play a crucial role in future work. There is evidence that even an exact second virial coefficient is not sufficient to predict a first order nematic–smecticA transition for l → ∞. Indeed, there are indications that also the third term of the present functional from Eq. (14), which only contributes to higher virial coefficients, is still captured insufficiently. First, the functional does not show the excellent Percus-Yevick direct correlation function for the hard sphere fluid. Second, the third virial coefficient is not reproduced exactly for all hard-body systems other than hard spheres. Finally, the isotropic–nematic transition in a system of hard thin platelets is best described after an additional semi-empirical modification. The presented dimensional crossover result and scaling arguments can be directly applied to obtain an

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-14

Wittmann, Marechal, and Mecke

expression in terms of higher rank tensors which might remedy some of these shortcomings. In the final part of Sec. III B, we have argued about the inconsistency between an exact second term and the representation of the third term with a finite number of tensors, which is related to the scaling of the orientational order parameter with the aspect ratio. An alternative representation for the third term would be desirable, as the direct implementation of the general kernel from Eq. (16) appears unrealistic. The ultimate test for an improved functional is the solid state. A functional which has the correct crossover properties to lower dimensions should describe the strongly confined particle assembly in a crystal very accurately. The high computational effort could be reduced by the restriction to parallel rods, which, using edFMT with an adequate expression for the third term, should be described in an equal way as with the functional for parallel HC.14, 15 Doing so, however, edFMT would still face the same problems as all fundamental measure functionals based on scaled particle theory that the fluid phases are not well described. Therefore, it is a difficult task for future work, to further improve the EOS and phase diagram of HSC, without spoiling the generality for arbitrary shapes. In conclusion, we can say that this work illustrates a substantial progress in understanding the relation between current approximations in edFMT and the virial series of anisotropic hard bodies. We provided the currently best expression of the functional for rods at higher densities and made the first steps towards a more general approach, which hopefully obviates empirical corrections. ACKNOWLEDGMENTS

The authors like to thank Roland Roth for many helpful suggestions. The system of hard cylinders was first considered in the Bachelor’s thesis by Sebastian Kaul. We also thank Hendrik Hansen-Goos for sharing the numerical values of ζ (S). Financial support by the DFG under Grant No. Me1361/12 as part of the Research Unit “Geometry and Physics of Spatial Random Systems” is gratefully acknowledged.

J. Chem. Phys. 141, 064103 (2014)

precisely, we set θ = 0 and obtain the LO term  2π   dφ (ν) (pa) ωC ◦ (θ = 0, φ) = O l 0 (A1) nν,C (z) = ρD C0 (z) 2π 0 of the weighted densities with the common integral  1 2   l C0 (z) := dzL ρC (z, zL ) = O l 1 ρ 0

(A2)

of all weighted densities and a one-dimensional density modulation ρ(z), combined to ρC (z, zL ) = ρ (zL L + z) + ρ (−zL L + z).

(A3)

The hemispherical contribution  D 1 (pa) nν,H (z) = dz ρ (ν) (z, zD ) 2 0 D H    2π   L dφ (ν) D × ωH◦ zD + , θ = 0, φ = O l −1 2π 2 2 0 (A4) to the weighted densities of the perfectly aligned system can (ν) be calculated with the weight functions ωH ◦ (z, θ, φ) from Eqs. (B3) to (B6) in Ref. 28. We simplified the integration using the combined density modulation     L L D D (ν) + ρ z − zD − ρH (z, zD ) = σν ρ z + zD + 2 2 2 2 (A5) with the substitution z → zD D2 + L2 and σ ν = −1 to calculate a vectorial weighted density and σ ν = 1 otherwise. The Taylor expansion of the full one-dimensional weight functions ω(ν) (z,  ) from Ref. 28 gives rise to the linearized orientation-dependent weighted densities     S ∂2 dφ (S )   (ν)   ρ(z )ω (z − z ,  ) . nν 2 (z) = 2 dz  2 3 ∂θ 2π θ=0 (A6)

APPENDIX A: LINEARIZATION OF THE FUNCTIONAL

The calculation  requires the assumption, that only small angles θ = O l −p occur between the particle axis and the (neˆ eˆz , for some exponent p > 0. The orientamatic) director n tional average becomes  π   2π 2 1 dφ θ dθ (A7) d  2π 0 0

It is not possible, to calculate inhomogeneous weighted densities with full orientational degrees of freedom analytically. The artificial system of strictly parallel bodies, however, can be evaluated easily. Inspired by the scaling behavior towards l → ∞, where only the leading order of the functional is relevant, we can further define weighted densities, which are linearized in the orientation dependence. As such, we find analytic expressions in Appendix B and present a free minimization of the resulting linearized orientational distribution. As a starting point we consider the one-dimensional weight functions of HSC. We see from Eq. (B21) in Ref. 28, that the LO term in the aspect ratio l is exclusively due to the cylindrical contributions ωC(ν)◦ (θ, φ), without explicit dependence on z, defined in Eqs. (B12)–(B16) therein. More

and from the Taylor expansion of cos 2 θ in Eq. (23), we find the NLO    3 S2 = θ dθ θ 2 g2 (θ ) = O l −2p (A8) 2   of the order parameter S = 1 − S2 + O l −4p . We need to introduce the orientational decoupling approximation to the scaled distribution g2 ( ) in Eq. (2), to accordingly define S2 in Eq. (A6), as the integration over dz and the Taylor expansion need to be carried out first. For now, we are particularly interested in the terms, which contribute to the limit l → ∞. As discussed in Sec. III A, the diverging LO term in the excess free en(TR) ergy per particle  that   0  vanishes if ζ = 2 and φ3 = φ3 , such ex /ρ = O l . Hence, the NLO contributions of O l −1 to

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-15

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

the weighted densities are still important at l → ∞. The expansion of the density   η 4η 8η (A9) ρ= = 3 − + O l −3 , v D πl 3D 3 π l 2 allows to extract the explicit LO and NLO terms from Eqs. (A1) to (A4). It remains to be seen, whether Eq. (A6) also provides a NLO contribution, as the exponent p in Eq. (A8) is yet undetermined. To this end, we neglect the hemispherical contribution and only consider the LO orientation-dependent term in (S )

nν 2 (z) = ρD

2S2 3



  dφ (ν) C2 (z, φ) + O S2 l −1 2π

(A10)

of the cylindrical weighted density with  ∂ 2 ωC(ν)◦ (θ, φ)  (ν) (ν) C2 (z, φ) = C0 (z)  + ωC ◦ (0, φ)Cθ (z) 2 ∂θ 2  θ=0 (A11) and   cos θ  2   l ∂2  Cθ (z) = dz ρ (z, z ) = O l 1 . (A12) L C L  2  2ρ ∂θ 0 θ=0

Note, that the perfectly aligned cylindrical weights from Eq. (A1) provide the zero order to the expansion in Eq. (A10). Now, we combine the linearized weighted densities (pa)

(S )

(pa)

n(lin) = nν,C + nν,H + nν 2 ν

(A13)

to the scaled functional (lin) ex and calculate the variation    δ ∂ex (z) 3 2 = dr ρ(z)(ln g2 (θ ) + 1) + θ θ =0 δg2 (θ ) ∂S2 2 (A14) to determine the orientational distribution function g2 (θ ) =

2α 2 e−α 

2 2

θ

 (A15) α2 1−e for a highly aligned system. We find for the intrinsic order parameter    1 d 3 ∂ex (z) (A16) dz = O l 1 , α 2 := d 0 2ρ ∂S2   as we deduce ∂ex /∂S2 = O l 0 from Eq. (A10). From Eq. (A8), we obtain   π2 2  − 4 α π2 2   31− 1+ 4 α e 3  + O l −2 , S2 = =  2  2 π 2 2α1 l − α2 α2 1 − e 4

A simple calculation yields the ideal gas term       id (z) = ρ(z) ln 2ρ(z) α12 − 2 + O l −2 consistently. Moreover, the explicit expression  (pa)  ex (z) ex (z) ˜ lim = lim + C(z)S2 l→∞ l→∞ ρ ρ

(A18)

(A19)

for the excess free energy becomes particularly simple, as there remains only a linear term in the correction S2 to the (pa) order parameter. The functional ex for the strictly parallel (pa) (pa) system is obtained with the weighted densities nν,C and nν,H and the prefactor C˜ can be easily determined in a similar way. With the corresponding result of Eq. (A16), we see from Eq. (A17), that S2 is exactly proportional to the inverse of the in˜ tegral over its prefactor C(z), leaving a trivial term. Hence, we find the overall free energy density F = F (pa) + 2N ln α1 + C

(A20)

which reminds of the expression found in Ref. 41. Only the orientational entropy 2Nln α, which arises from Eq. (A18)  with the overall particle number N = drρ(z), explicitly depends on the order parameter. Other contributions of the ideal gas contribute to the overall free energy F (pa) of perfectly aligned bodies or are absorbed into the constant C, which is not important for the phase behavior. At finite aspect ratio, the functional does not take the simple form of Eq. (A20), due to contributions, which are quadratic, cubic or reciprocal in S2 . Recall, that we also need to consider the full linear dependence from Eq. (A6) on this order parameter in Eq. (A13). In a good approximation, we may set the limit of the integration over θ in Eq. (A7) to infinity, as the result for α in Eq. (A16) ensures a quick de2 2 cay of the orientational distribution g˜ 2 (θ ) = 2α 2 e−α θ at high enough aspect ratio. Doing so, only the simple leading order S2 = 3/(2α 2 ) from Eq. (A17) remains.

2 − π4

(A17) where we only use the LO α 2 = α12 l in the final step. This result is in accordance  1  with Eq. (26) and gives rise to p = 1/2 and thus θ = O l 2 , as derived in Eq. (32). Now, we see that   (S ) the weighted densities nν 2 = O l −1 from Eq. (A10) provide an additional NLO term to the perfectly aligned weighted (pa) (pa) densities nν,C and nν,H , while all other orientation-dependent terms scale out in the limit l → ∞.

APPENDIX B: WEIGHTED DENSITIES OF ALIGNED SPHEROCYLINDERS

The weighted densities n(lin) ν , defined in Eq. (A13) of Appendix A, are sufficient to appropriately describe perfectly aligned systems at any aspect ratio l and orientation dependence in the limit l → ∞. Here, we provide the explicit expression using the cosine modulation from Eq. (47). The results π (pa) n3,C = ρ D 3 C0 (z), 4 ← (pa) (pa)  n = ρπ D 2 C (z), (B1) n =2 → 2,C

(pa)

n1,C

2,C 11

0

 →(pa)   →(pa)  1 =2 ← n 1,C 11 = − ← n 1,C 33 = ρ DC0 (z) 4

for the LO, defined in Eq. (A1), can be written with the common cylindrical integral C0 (z) =

d Lπ L π + (z) sin = l + (z)L sin , D Dπ d dL (B2)

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-16

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

where dL := d/L, L := ldL /π , and (z) = λ cos 2πz . We see d  0 (pa) that nν,C = O l . The hemispherical contributions   1 (z)  π = + L − 23L sin ρ D3π 6 4 dL (pa)

n3,H

and the hemispherical integrals n    n 2g+1 (2g + 1)! g L Hn (z) = (z) g



over the density modulation, where

π (1 + l) π (1 + l) − 22L cos dL l dL l   (pa) n2,H π (1 + l) π = 1 + (z) L sin − L sin ρ D2π dL l dL =2

(pa) n0,H

= ρD ρ    → (− n 2,H )3 π (1 + l) π d 2L sin = z− − 2L sin ρ D2π 4 dL l dL   (pa) − →  n 1,H π (1 + l) 3 =2 − cos dL l ρD ←  (pa)  → n 2,H 11 1 (z) π = + (−L − 23L ) sin 2 ρD π 3 2 dL  π (1 + l) π (1 + l) 3 2 − 2L cos + 2L sin dL l dL l ←  (pa)  → n 2,H 33 1 π (1 + l) = + (z) (L − 23L ) sin ρ D2π 3 dL l  π π (1 + l) + 23L sin (B3) + 22L cos dL dL l   (pa) from Eq. (A4) are nν,H = O l −1 . Substituting the density with Eq. (A9), we find the explicit leading terms (pa)

  4η Cθ (z) + O l −2 , 2 3l   8η − 2 Cθ (z) + O l −2 , l D   4η − 2 C (z) + O l −2 , l π D2 θ   8η (B4) − 2 C (z) + O l −2 , l π D3 θ     d π 2η z− sin + O l −2 , lD 4 dL     π η d sin z− + O l −2 , lπ D 2 4 dL     8η (pa) ← → n 2,H = − 2 Cθ (z) + O l −2 , 33 3l D

n3,H = − (pa)

n2,H = (pa)

n1,H = (pa)

n0,H = − (pa)  → n 2,H 3 = − (pa)  → n 1,H 3 = ← (pa)  → n 2,H 11 =

where (z − d4 ) = λ sin 2πz and d   π l . 1 + (z) cos Cθ (z) = − 2 dL

(B5)

[2n] As calculated in Appendix C the product [n[2n] 1 n2 ] of perfectly aligned tensorial weighted densities of rank 2n in Eq. (49) consists of the cylindrical integral C0 from Eq. (B2)

π −2h π (1 + l)  (−1)h L − sin dL l h=0 (2h)! dL g

+ 23L sin

(pa) n1,H

(B6)

g=0

g = sin

−(2h+1) π (1 + l)  − cos (−1)h L . dL l h=0 (2h + 1)! g−1

(B7)

The exact hemispherical contribution of all tensors in n12 from Eq. (51) reads    π π (1 + l) H = (z) 22L cos − cos dL dL l     π (1 + l) π −L sin = O l −1 . (B8) + sin dL l dL Interestingly, the contribution ρH is of lower order than NLO, in contrast to the hemispherical weight functions from Eq. (B5). This can be generally shown for Hn in Eq. (B6). Defining ∞  −2h ¯ g = sin π (1 + l) (−1)h L  dL l h=g+1 (2h)!

− cos

∞ π (1 + l)  −(2h+1) , (−1)h L dL l h=g+1 (2h + 1)!

(B9)

we find

n     Hn (z)  n  ¯ g 2g+1 (2g + 1)! ∞ −  = L g (z) g=0    π (1 + l) g (−1) = O l −1 (B10) + cos dL l

as  ∞ includes series representations of the trigonometric functions and vanishes after application of an addition theorem. The last expression, obtained from subtracting the term with h = g, becomes zero after summation over g. The orientation-dependent cylindrical NLO terms   4ηS2 (S ) (C0 (z) + 2Cθ (z)) + O l −2 n2 2 = 3lD 4 (S2 ) (S ) n3 = 4Dn1 2 = D ←    4ηS (S ) → 2 Cθ (z) + O l −2 n 2 2 11 = 3lD (B11) ←   −2  4ηS (S2 ) → 2 C (z) + O l n 2 33 = 3lD 0 ←   ηS2 (S )  → n 1 2 11 = − (C0 (z) − Cθ (z)) + O l −2 2 3lπ D 1  →(S2 )  n 1 33 , = − ← 2   (S ) defined in Eq. (A10), also are nν 2 = O l −1  . Note that for the nematic phase we have C0 + 2Cθ λ=0 = 0. All

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-17

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

contributions to the scaled weighted densities vanish, if they are not listed explicitly. APPENDIX C: TENSOR SERIES FOR PARALLEL HARD SPHEROCYLINDERS

In this appendix, we explicitly calculate the product    [r] [r] (12)   n1 n2 = dz dz ρ (z + z ) ρ (z + z ) ω[r] (z , z )  D3  (I )  I ) (I I )  ρ 2 16π C0 n C0 + Hn(I ) + (I n Hn (C1) = 0 if r = 2n of two tensorial weighted densities of arbitrary rank r in Eq. (49), using cylindrical coordinates r = (z, , ϕ). The combined weight function    (1) (12) (2)     ω[r] (z, z ) =  dd (r) ω[r] (r ) dϕdϕ  ω[r] (C2) of parallel HSC implies the contraction of tensorial weight functions from Eq. (6). The trivial integrals over  and  yield a prefactor D2 /4. The z integrands for the aligned cylindrical parts only include the density modulation, as the unit vectors ⎞ ⎞ ⎛ ⎛ ⎛ ⎞ cos ϕ − sin ϕ 0 ⎟ I ⎜ ⎟ II ⎜ ⎟ ⎜ (C3) n = ⎝ sin ϕ ⎠ v = ⎝ cos ϕ ⎠ v = ⎝ 0 ⎠ 0

0

1

of an aligned cylinder are constants. The respective results are proportional to C0 as defined in Eq. (A2). The weight function (1) ω[r] only has such a contribution with the deviatoric curvature κ(z) = 1/D of a cylinder and κ(z) = 0 for a sphere. Hence, the only non-trivial integration over z includes the hemispherical normal vector  ⎛ ⎞ 2 cos ϕ 1 − zD ⎜  ⎟ (C4) n=⎜ 2 ⎟ ⎠ ⎝ sin ϕ 1 − zD zD (2) in ω[2n] , with the substitution z → zD D2 + L2 introduced in (1) Appendix A for the upper hemisphere. As all elements of ω[2n] which include at least one third component of n or vI vanish, there are only two types of integrals over z in Eq. (C1). In analogy to Eq. (A4) we can calculate  1   1 (I ) 2 n dzD ρH (z, zD ) 1 − zD (C5) Hn (z) = 2ρ 0

and Hn(I I ) (z) =

1 2ρ

 0

1

  2 (n−1) 2 dzD ρH (z, zD ) 1 − zD zD

(C6)

from the normal vector in Eq. (C4) and the combined density from Eq. (A5). The first integral type only considers the first two vectorial components and is related to the principal cur(1) ) vature direction vI in ω[2n] . Thus the factor (I n in Eq. (C1) collects all contributing integrals over ϕ and ϕ  for the components α 1 . . . α r in Eq. (C2), which have the form  2π dϕ sinh ϕ cosk ϕ , h + k = r . (C7) 0

The hemispherical integral Hn(I I ) from Eq. (C6) corresponds to α 1 . . . α r − 2 ∈ {1, 2} and α r − 1 = α r = 3, i.e., the principal curvature direction vII . In this case, there is no cylindri(2) I) proportional to (I in Eq. (C1). cal contribution from ω[2n] n This term includes the sum over integrals of the same type as Eq. (C7), but with h + k = r − 2. For mixed combinations including both principal curvature directions or odd tensor ranks, we find at least one odd value of h or k in all relevant tensor components. In conclusion, we see that only even tensor ranks r = 2n with the prefactors   I) 1 (I 3  n − n (I ) 2 = 2π 2 (C8) n = n = 2n (n + 1) do not vanish, as the integrals over odd powers of sin ϕ or cos ϕ in Eq. (C7) are zero. Combination of Eqs. (C5) and (C6) and explicit calculation of the integrals for the density profile given by Eq. (47) results in Eq. (49) with the combined contribution Hn := Hn(I ) + 2nHn(I I ) from Eq. (B6). However, the calculation of the exact contribution to the second virial coefficient from Eq. (C1) is not possible, as the sum over all tensorial weighted densities cannot be evaluated. We can calculate the contribution  1 1 dzD ρH (z, zD ) H= 2ρ 0 ×

∞      n  2 n 2 2 (n−1) 1 − zD 1 − zD − 2nzD 2 4π n=1

1 = 2ρ



0

1

dzD ρH (z, zD ) (1 − 2|zD |)

(C9)

to Eq. (51) without the mathematically questionable detour of processing individual tensorial weighted densities in Eqs. (C1) and (C2), when we interchange the integrations with the subsequent summation over n in the spirit of Eq. (50). This is because the series  of weighted densities involves an evaluation of the sum ∞ n=1 n Hn after the spatial integrals, although the expansion of Eq. (10) is performed within the integral. With the density modulation in Eq. (47) we obtain H as stated in Eq. (B8). For a homogeneous density, i.e., the nematic phase, a quick calculation shows that  [2n] [2n] L2 D  , = ρ2 n1 n2 16π n

(C10)

as all Hn and the exact H vanish. In this case, we can eas [2n] ily verify the equivalence n [n[2n] 1 n2 ] = n12 between the series of tensors in Eq. (49) and the approach from Eq. (51).

APPENDIX D: EXCLUDED VOLUME OF TWO SPHEROCYLINDERS

An argumentation analogous to the one in Appendix C yields the tensorial contribution   (1) 1 (2) gn (γ ) := − 2 (r, 0)ω[2n],R (r , γ ) drdr ωR,[2n] L D (D1)

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-18

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014)

of rank 2n to the mutual excluded volume   ∞ π  4 2 vFMT (γ ) = 2L D + gn (γ ) + 2π LD 2 + π D 3 4 3 n=1 (D2) of HSC. The subscript R denotes orientation dependence of the unit vectors within the weight function as introduced in Ref. 28. As there is no density modulation in the spatial integrands, the hemispherical contributions vanish as for the nematic phase in the parallel system. The cylindrical weight functions, however, are now orientation-dependent. More specifically it is sufficient to account for the intermolecular angle γ by respectively choosing the arguments θ = 0 (ν) (r, θ ) in and θ  = γ of the rotated weight functions ωR,[2n] Eq. (D1). As described in Ref. 28, we need to use the rotated normal vector ⎛ ⎞ cos ϕ cos γ ⎜ ⎟ sin ϕ n(ϕ) = ⎝ (D3) ⎠ − cos ϕ sin γ

n (γ ) = −16π

Πnk cos γ =: −16πgn (γ ) 2k

(D4)

k=0

with the coefficients (2n − 2k)!(2k + 1)! π (2n)! (D5) Πnk = 2n+1 2 4 n! (2n − 1) (n − k)!2 k!2 (2n − 2k − 1) depends on the intermolecular angle γ , such that n (0) = n from Eq. (C8). The result for gn (γ ) is closely related to Eq. (C10). The overall tensorial contribution ∞ 

gn (γ ) = −Π00 +

∞  ∞ 

Πnk (cos2 γ )k

k=0 n=k

n=1



=−

sin γ −

π π  (2k − 2)!(cos2 γ )k − = sin γ − 2k−1 (k − 1)!k! 4 2 4 k=0 (D6)

can be calculated after interchanging the sums over n and k. This result shows that the exactly known expression for the mutual excluded volume, stated in Eq. (37), is reproduced with edFMT via Eq. (D2). However, it is not possible to evaluate sum over k in Eq. (D4) and make an analytic statement about convergence of Eq. (D6).

∞  π d2n P2n (cos γ ) =− 4 n=1

(E1)

in Eq. (D6) can be expanded in terms of even Legendre Polynomials P2n (cos γ ) with the coefficients58 d2n = π

(4n + 1)(2n − 3)!!(2n − 1)!! . 22n+2 n!(n + 1)!

(E2)

We introduce the prefactors ζm[n] to rank 2m tensors, such that the respective edFMT contributions −

to include the orientation dependence, where we set φ = 0 without loss of generality. Rotation by θ = 0 maintains the unit vectors from Eq. (C3) for the first weight function. Hence, we only need to consider the two cases, described in Appendix C, which now have orientation-dependent prefac(I /I I ) (γ ) from integrals similar to Eq. (C7). The comtors n bined expression n 

tional, evaluated in Eq. (D6), converges slowly and truncation after rank six tensors does not clearly improve upon the ζ corrected rank two tensors from Eq. (11), as recognized in Ref. 27. Now, we derive the generalized semi-empirical cor(ζ [n ]) rection for additional higher rank tensors in φ2 t from Eq. (12), which is suited for a sophisticated description of the nematic phase. The result

n 

ζm[n] gm (γ ) = d2n P2n (cos γ )

(E3)

m=1

to the excluded volume add up to the term of order n of the expansion in Eq. (E1). Inserting the general expressions for gn from Eq. (D4), we obtain the explicit expression for ζm[n] in Eq. (13). The first values read 5 3 21 , ζ [2] = − , ζ2[2] = , 4 1 8 8 13 117 429 , ζ [3] = − , ζ3[3] = . ζ1[3] = 64 2 32 64 The generalized ζ correction ζ1[1] =

∞ 

n[2m] n[2m] 1 2



m=1

nt n  

ζm[n] n[2m] n[2m] 1 2

(E4)

(E5)

n=1 m=1

may be readily applied to the nematic weighted densities. This new approximation improves on pure truncation of the left-hand side and both expressions are exact in the limit n nt ζm[n] of nt → ∞. Note that the overall prefactor ζmt := n=m each term effectively depends on the maximal tensor rank n 2nt and ζmt → 1 for nt → ∞. For parallel HSC we (ζ ) can verify n = 16π d2n in Eq. (35) by definition of Eq. (E3) for γ = 0. The construction of the generalized ζ correction from Eq. (E3) further ensures that each sum n 

ζm[n] n[2m] n[2m] = d2n P2n 1 2

(E6)

m=1

assigned to tensors of rank 2n constitutes one order parameter  1 P2n = dxP2n (x)g [nt ] (x) (E7) 0

APPENDIX E: CORRECTION OF HIGHER ORDER TERMS

In Appendix D, we calculated the general contribution gn (γ ) of rank 2n tensors to the mutual excluded volume of HSC. The infinite series in the second term φ2[∞] of the func-

to the excess free energy of a homogeneous nematic. A straight forward minimization yields the orientational distribution function  n  t  [nt ] 2 α[2n] P2n (cos θ ) (E8) g (cos θ ) = C exp n=1

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

064103-19

Wittmann, Marechal, and Mecke

J. Chem. Phys. 141, 064103 (2014) 28 R.

of order nt with the intrinsic order parameters 2 α[2n] := −

2 d2n P2n ∂ed = + A3 ({P2m }) , ρ ∂P2n 1−η

(E9)

where n, m ≤ nt . The second term A3 , which couples the equations arises from the third term φ 3 of the functional and depends on its explicit form. The nematic equilibrium state is obtained after the self-consistent calculation of the normalization factor C in Eq. (E8) and solution of nt equations for the nematic order parameters from Eq. (E8). In the Onsager limit, our generalized ζ correction exactly reproduces the results from Ref. 58 by construction. This expansion is reported to converge at nt = 7 for isotropic–nematic coexistence, which our calculation confirms for finite aspect ratios. 1 R.

Evans, Adv. Phys. 28, 143 (1979). Rosenfeld, Phys. Rev. Lett. 63, 980 (1989). 3 J. K. Percus and G. J. Yevick, Phys. Rev. 110, 1 (1958). 4 H. Reiss, H. L. Frisch, E. Helfand, and J. L. Lebowitz, J. Chem. Phys. 32, 119 (1960). 5 Y. Rosenfeld, M. Schmidt, H. Löwen, and P. Tarazona, J. Phys.: Condens. Matter 8, L577 (1996). 6 Y. Rosenfeld, M. Schmidt, H. Löwen, and P. Tarazona, Phys. Rev. E 55, 4245 (1997). 7 P. Tarazona and Y. Rosenfeld, Phys. Rev. E 55, R4873 (1997). 8 P. Tarazona, Phys. Rev. Lett. 84, 694 (2000). 9 M. Marechal, H. H. Goetzke, A. Härtel, and H. Löwen, J. Chem. Phys. 135, 234510 (2011). 10 J. A. Cuesta, Phys. Rev. Lett. 76, 3742 (1996). 11 J. A. Cuesta and Y. Martínez-Ratón, Phys. Rev. Lett. 78, 3681 (1997). 12 J. A. Cuesta and Y. Martínez-Ratón, J. Chem. Phys. 107, 6379 (1997). 13 Y. Martínez-Ratón and J. A. Cuesta, J. Chem. Phys. 111, 317 (1999). 14 Y. Martínez-Ratón, J. A. Capitán, and J. A. Cuesta, Phys. Rev. E 77, 051205 (2008). 15 J. A. Capitán, Y. Martínez-Ratón, and J. A. Cuesta, J. Chem. Phys. 128, 194901 (2008). 16 Y. Martínez-Ratón and J. A. Cuesta, Phys. Rev. Lett. 89, 185701 (2002). 17 Y. Martínez-Ratón and J. A. Cuesta, J. Chem. Phys. 118, 10164 (2003). 18 Y. Martínez-Ratón, Phys. Rev. E 69, 061712 (2004). 19 M. Schmidt, Phys. Rev. E 63, 050201(R) (2001). 20 J. M. Brader, A. Esztermann, and M. Schmidt, Phys. Rev. E 66, 031401 (2002). 21 A. Esztermann, H. Reich, and M. Schmidt, Phys. Rev. E 73, 011409 (2006). 22 L. Onsager, Ann. NY Acad. Sci. 51, 627 (1949). 23 G. Cinacchi and F. Schmid, J. Phys.: Condens. Matter 14, 12223 (2002). 24 Y. Rosenfeld, Phys. Rev. E 50, R3318 (1994). 25 Y. Rosenfeld, Mol. Phys. 86, 637 (1995). 26 H. Hansen-Goos and K. Mecke, Phys. Rev. Lett. 102, 018302 (2009). 27 H. Hansen-Goos and K. Mecke, J. Phys.: Condens. Matter 22, 364107 (2010). 2 Y.

Wittmann and K. Mecke, J. Chem. Phys. 140, 104703 (2014). Marechal and H. Löwen, Phys. Rev. Lett. 110, 137801 (2013). 30 M. Marechal, U. Zimmermann, and H. Löwen, J. Chem. Phys. 136, 144506 (2012). 31 S. C. McGrother, D. C. Williamson, and G. Jackson, J. Chem. Phys. 104, 6755 (1996). 32 P. Bolhuis and D. Frenkel, J. Chem. Phys. 106, 666 (1997). 33 J. M. Polson and D. Frenkel, Phys. Rev. E 56, R6260 (1997). 34 A. Poniewierski, Phys. Rev. A 45, 5605 (1992). 35 A. Poniewierski and R. Hołyst, Phys. Rev. Lett. 61, 2461 (1988). 36 R. Hołyst and A. Poniewierski, Phys. Rev. A 39, 2742 (1989). 37 A. Poniewierski and R. Hołyst, Phys. Rev. A 41, 6871 (1990). 38 A. Poniewierski and T. J. Sluckin, Phys. Rev. A 43, 6837 (1991). 39 A. M. Somoza and P. Tarazona, Phys. Rev. Lett. 61, 2566 (1988). 40 A. M. Somoza and P. Tarazona, J. Chem. Phys. 91, 517 (1989). 41 A. M. Somoza and P. Tarazona, Phys. Rev. A 41, 965 (1990). 42 E. Velasco, L. Mederos, and D. E. Sullivan, Phys. Rev. E 62, 3708 (2000). 43 J. D. Parsons, Phys. Rev. A 19, 1225 (1979). 44 S.-D. Lee, J. Chem. Phys. 87, 4972 (1987). 45 M. Esposito and G. T. Evans, Mol. Phys. 83, 835 (1994). 46 Y. Martínez-Ratón, E. Velasco, and L. Mederos, J. Chem. Phys. 123, 104906 (2005). 47 H. Graf and H. Löwen, J. Phys.: Condens. Matter 11, 1435 (1999). 48 H. Graf and H. Löwen, Phys. Rev. E 59, 1932 (1999). 49 A. Stroobants, H. N. W. Lekkerkerker, and D. Frenkel, Phys. Rev. Lett. 57, 1452 (1986). 50 A. Stroobants, H. N. W. Lekkerkerker, and D. Frenkel, Phys. Rev. A 36, 2929 (1987). 51 J. A. C. Veerman and D. Frenkel, Phys. Rev. A 43, 4334 (1991). 52 B. Mulder, Phys. Rev. A 35, 3095 (1987). 53 R. van Roij, P. Bolhuis, B. Mulder, and D. Frenkel, Phys. Rev. E 52, R1277 (1995). 54 N. F. Carnahan and K. E. Starling, J. Chem. Phys. 51, 635 (1969). 55 R. Roth, R. Evans, A. Lang, and G. Kahl, J. Phys.: Condens. Matter 14, 12063 (2002). 56 H. Hansen-Goos and R. Roth, J. Phys.: Condens. Matter 18, 8413 (2006). 57 G. J. Vroege and H. N. W. Lekkerkerker, Rep. Prog. Phys. 55, 1241 (1992). 58 H. N. W. Lekkerkerker, P. Coulon, R. Van Der Haegen, and R. Deblieck, J. Chem. Phys. 80, 3427 (1984). 59 R. L. C. Vink, S. Wolfsheimer, and T. Schilling, J. Chem. Phys. 123, 074901 (2005). 60 S. Wolfsheimer, C. Tanase, K. Shundyak, R. von Roij, and T. Schilling, Phys. Rev. E 73, 061703 (2006). 61 F. Schmid, G. Germano, S. Wolfsheimer, and T. Schilling, Macromol. Symp. 252, 110 (2007). 62 E. Velasco, L. Mederos, and D. E. Sullivan, Phys. Rev. E 66, 021708 (2002). 63 K. Shundyak and R. van Roij, J. Phys.: Condens. Matter 13, 4789 (2001). 64 D. Frenkel and R. Eppenga, Phys. Rev. Lett. 49, 1089 (1982). 65 M. Marechal, A. Cuetos, B. Martínez-Haya, and M. Dijkstra, J. Chem. Phys. 134, 094501 (2011). 66 R. Roth, J. Phys. Condens. Matter 22, 063102 (2010). 67 M. S. Wertheim, Mol. Phys. 83, 519 (1994). 29 M.

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 128.114.34.22 On: Mon, 01 Dec 2014 15:02:24

Fundamental measure theory for smectic phases: scaling behavior and higher order terms.

The recent extension of Rosenfeld's fundamental measure theory to anisotropic hard particles predicts nematic order of rod-like particles. Our analyti...
893KB Sizes 1 Downloads 7 Views