Analytical Virial Coefficients and New Equations of State of Hard Ellipsoid Fluids
Maryam Hashemi, Abolghasem Avazpour*, Shaker Hajati
Department of Physics, Yasouj University, Yasouj, Iran
To cite this article:
Maryam Hashemi, Abolghasem Avazpour, Shaker Hajati. Analytical Virial Coefficients and New Equations of State of Hard Ellipsoid Fluids. Chemical and Biomolecular Engineering. Vol. 2, No. 1, 2017, pp. 27-40. doi: 10.11648/j.cbe.20170201.15
Received: January 2, 2017; Accepted: January 19, 2017; Published: February 22, 2017
Abstract: The complexity of calculations for high order virial coefficients of ellipsoids makes it difficult to obtain accurate analytical high order coefficients and equation of state for such systems. Using analytical method, the virial coefficients up to third order are calculated. For higher ones, the numerical values were taken from publications of other researchers, based on Monte Carlo integration method. By fitting the available numerical virial coefficients, sixth to eighth order, the two shape parameter analytical expressions of the hard convex molecules are obtained. Using these available data, up to eighth order, we have obtained the approximate one shape parameter analytical expressions of the hard prolate and oblate ellipsoid molecules. The fitted virial coefficients are in agreement with the simulation results. Moreover, the approximate analytical expressions for the equation of state of isotropic hard ellipsoid fluids are proposed. The proposed equations of state are in good agreement with the simulations up to medium elongations. In addition, our equations show a better agreement comparing to other works. Also, the newest equation is used for both prolate and oblate ellipsoid fluids and is convenient for elongations, k<10.0.
Keywords: Virial Coefficients, Equation of State, Hard Ellipsoid, Isotropic Fluid, Prolate Molecule, Oblate Molecule
Hard sphere is the simplest system in the statistical mechanics and is applied as a reference system for perturbation methods. This system is defined with interaction potentials that is considered to be repulsive forces between molecules . The simplicity of this model caused thermodynamic properties to be derived, using theoretical or computer simulation methods. The equation of state (EOS) of a system is its most important relationship between thermodynamic quantities such as pressure, volume and temperature, as it allows one to calculate most of its thermodynamic properties . Accurately equation of state is a thermodynamics equation that describes the status of the material under a set of special physical conditions. This equation that is derived by using mathematical relation between two or some of thermodynamic quantities is related to material type. Most of the prominent use of equation of state is the prediction of phase of gases and liquids. The simplest equation of state is ideal gas one is reasonably accurate at low pressures and intermediate temperatures. Moreover, many accurate equations of state have been developed .
The first convenient equation of state for liquid and gas phases was proposed by Van der Waals in 1873 [3,4]. Van der Waals deduced his equation with intuition. In 1881 Clausius modified Van der Waals’s equation. Thiesen suggested unlimited power series for non-ideal behavior of real fluids in 1885. These series can be expressed in terms of powers of density. The coefficients of the series were named virial coefficients . The importance of these coefficients lies in the fact that they are related directly to the interaction between
molecular clusters. Since virial series converged slowly, researchers used approximations such as Levin and Pade to accelerate convergence . About 80 equations of state for hard sphere systems have been proposed. Some of these equations were reported in Refs. [7-16].
Due to molecules anisotropy of most real liquids, liquid crystals and molecular fluids with convex particles have been studied. System of hard particles such as hard disk, hard sphere, hard rods and hard ellipsoids are appropriate models for simulation and theoretical studies . Hard spheres are models for real atomic liquids. Sphero- cylinders and uniaxial ellipsoids of revolution models are suitable for real fluids including diatomic or polyatomic linear molecules. Biaxial ellipsoids are models for nonlinear polyatomic molecules . Since the best approximation for molecules such as 4-Cyano-4'-pentylbiphenyl (5CB) and 4-Octyl-4’-Cyanobiphenyl (8CB) is ellipsoidal model, it is interested to study ellipsoidal molecular fluids.
Various equations of state have been suggested for hard ellipsoidal molecular liquids . The compressibility factor has been expressed by using packing fraction and shape parameter. Such equations are Boublik , Nezbeda  and Song -Mason  equations of state. The results of these equations for convex shapes, especially ellipsoids, show that using only one parameter for describing molecular anisotropy cannot be sufficient. Naumann and coworkers  considered a second shape parameter that could distinguish the prolate and oblate molecules. Other methods for obtaining the equation of state are based on resuming virial series. Barboy and Gelbart  proposed the compressibility factor as a polynomial function with modified variable that depends on volume fraction. Wojcik and Gubbins  also proposed a similar expression. Parsons  proposed a new equation of state based on the second virial coefficient of ellipsoids and the combination of spherical Carnahan- Starling equation.
A powerful systematic method for calculating bulk properties is via a virial expansion. In general, the virial coefficients depend on the temperature and the intermolecular potential. Using analytical methods, low order virial coefficients are calculated for simple intermolecular potentials. Therefore, virial coefficients up to fourth order were calculated for hard spheres [26,27]. The coefficients up to seventh order were calculated for hard parallel squares and cubes analytically . There are explicit expressions for the second virial coefficient of convex bodies . High order virial coefficients must be calculated numerically for the simplest interaction potentials.
The Monte Carlo method can also be readily applied to spherical hard bodies, and to date of the first eight coefficients have been reported for spheroids [30-33], prolate sphero-cylinders [32,34] and truncated spheres [32,35]. First six virial coefficients have been calculated for the hard Gaussian overlap (HGO) model  and the first five coefficients were calculated for oblate sphero-cylinders  and hard diatomics . First seven virial coefficients were calculated for a mixture of hard spheres and hard discs [32,39]. Rigby  has determined the first five virials of hard ellipsoids of revolution. According to our knowledge, there was not any analytical expression for sixth to eight virial coefficients of hard convex bodies in literature. Our motivation is to find analytical expressions for these coefficients and study their effects on the equation of state of hard ellipsoid fluid. Also, several equations of state proposed for hard convex bodies are tested.
In Section 2, the second and third-order virial coefficients of non-spherical molecules are introduced and related to Mayer function. Approximate analytical expressions for sixth to eighth virial coefficients in isotropic phase of the hard ellipsoids were obtained by using the shape parameters of convex molecule. In Section 3, approximate analytical expressions in the isotropic phase for the fourth up to eighth virial coefficients of hard prolate and oblate ellipsoids, are obtained versus length to breadth ratio of molecules and results are compared with the Monte Carlo data.
In Section 4, several equations of state proposed for hard convex bodies are introduced. The effect of virial coefficients on EOS, are investigated and its results are presented. Two new equations of state of hard ellipsoid liquid are proposed by using deduced virial coefficients up to eighth. The results are compared with the simulation data. The results, discussions and conclusions are presented in Sections 5 and 6.
2. Virial Coefficients of Molecular Liquids
The energy virial series for one component system is written as :
where , is the grand potential, contains the non-translational contributions to the ideal gas Helmholtz energy (i.e. arising from rotations, vibrations and electronic excitations) and is the chemical potential. Also, , and are the external potential, the de Broglie thermal wavelength and the one particle density(dependent in general both on position and orientation), respectively. is shorthand for the position, , and orientation, , of particle 1. The coefficient is given by
where is the nth virial coefficient for n particles with fixed positions and orientations. The second and third virial coefficients for spherical system are obtained by using statistical mechanics and linear algebra .
The interaction potential,
is dependent on the separation of two particles, , and Mayer function , is represented as:
In molecular liquids, the potential energy between molecules (1) and (2), , is function of the relative distances between molecular centers, , and directions of the axes of the molecules, . If the molecule is linear, , where , are the usual polar angles; if it is non-linear, , where , , are the Euler angles. Hard Gaussian overlap model are presented as a simple model for the theoretical study of liquid with a spherical molecules and uses in computer simulations. Berne and Pechukas  suggested this model for investigating interactions between ellipsoidal molecules, where strength of repulsive interaction between the molecules and is dependent on direction of molecules. The HGO potential for the two molecules and is
where molecular anisotropy is equal to:
is the range parameter with . For molecules with length to breath ratio of and unit length , anisotropy parameter is as follows:
Rickayzen  modified the Berne- Pechukas’s expression and suggested the distance parameter for two hard ellipsoids.
and are small and large semi-axes of ellipsoid.
Since the intermolecular potential depends on direction and relative distance of molecules, virial coefficients of ellipsoids are related to direction of molecules. For convex molecules, the second and third virial coefficients are written as follows:
where is the orientation distribution function. For isotropic liquids, the orientational distribution function, , and for the nematic liquids, , which minimizes the Onsager density functional, can be obtained from the self-consistence equation 
where is the normalization constant.
Using the second virial coefficient and hard Gaussian overlap model with closest approach of Rickayzen, is as follows:
The third virial coefficient, , depends on three particle interactions.
The second and third virial coefficients of hard ellipsoids of elongations; were calculated analytically . High order virial coefficients have more complicated integrals than the second and third coefficients. Isihara and coworker  showed that the second virial coefficient of hard convex bodies is changed significantly by the shape of particle. They showed that for any convex body, the average radius, surface area and volume are needed. Using calculations of Ref. , these parameters are obtained as
and , (19)
In Eqs. (16)- (18), are semi-axes of ellipsoid and the variable as a function of and is
After occurring these changes in variable, the available integrals are calculated from Gradshteyn and Ryzhik’s table of integrals and the results are expressed as :
and are the first and second kind of elliptic integrals given by
are defined as follow:
where is the single-particle density. By expanding Eq. (29) in terms of density, second virial coefficient is obtained as follows:
This result is consistent with the values obtained by Isihara and coworker . For hard convex bodies, this coefficient can be rewritten as follows:
where , the nonsphericity of the molecule is:
, and are the mean radius of curvature, surface and volume of hard ellipsoid respectively.
Because of the complex integrals in the higher virial coefficients of ellipsoids, these coefficients are calculated numerically. One of the most suitable numerical methods is application of Monte Carlo integration method . Computer simulations of hard ellipsoids were pioneered by Viellard – Baron  in two dimensions, and by Frenkel and Mulder  in three dimensions. Allen  performed computer simulations of biaxial hard ellipsoids. Researchers have calculated first eight coefficients of ellipsoids by using Monte Carlo method [31-33, 53-55].
When two of the main semi- axes of the ellipsoids have the same length, the molecule has uniaxial symmetry. Uniaxial hard ellipsoids of revolution are commonly denoted as prolate spheroids, when the third semi-axes is larger than the other two, and oblate spheroids, when it is smaller. If each of the three semi-axis has a different size, molecule will have biaxial symmetry. Prolate spheroids can be written as and oblate spheroids as, . Length to width ratio, , for the prolate spheroid is and for the oblate one is .
Carlos Vega  showed that the accuracy of higher virial coefficients of hard ellipsoids depends on two non- sphericity parameters. He had fitted successfully the data of third up to fifth virial coefficients by using two shape parameters. For hard convex bodies, he considered Minkowsky inequalities  given by
Using these inequalities, can be defined two following parameters:
By comparing these two definitions and definition of it is easy to show that
For a sphere, and are unity and for any other convex body and are larger than unity. The parameter was also used by Naumann and Leland  before in their research for a general equation of state of hard convex bodies. The value of for a prolate spheroid is identical to that of an oblate spheroid. However, the prolate spheroid presents a large value of than the oblate one . Vega analytical expressions have expressed in terms of two parameters and . These two new parameters are defined as follows:
After Vega's procedure, the expressions were fitted on the sixth to eighth virial coefficients of hard ellipsoid on the isotropic phase that was obtained by the Monte Carlo integration method of Refs. [32,33,53]. Our fitted expressions are as follows:
Reduced virial coefficient is defined as , where is the order of virial coefficient and is the volume of ellipsoidal molecule. In Tables 1- 3 sixth to eighth data obtained from the fitted expressions (41)- (43) are compared to numerical results of Monte Carlo method [32, 33, 53].
3. Approximate Expressions for the Virial Coefficients in Terms of Anisotropy of the Ellipsoidal Molecule
For hard uniaxial ellipsoid molecules, simple expressions require to express the virial coefficients. Thus, the length to width ratio of the molecules,, was used and the fitted expressions for the fourth to eighth virial coefficients are obtained. By using the numerically obtained data from Ref. , the one shape parameter fitted expressions were obtained:
To obtain the approximate analytical expressions for the sixth and seventh virial coefficients in isotropic phase, the data taken from Ref.  were used. These data were obtained by the Monte Carlo integration method. The following fitted expressions have been obtained.
|B6* (Fitted)||B6* (Simulation)|
|B7* (Fitted)||B7* (Simulation)|
Also, in Figures. 1-8, the to results are compared with the Monte Carlo data. All the virial coefficients data are presented in Refs. [32, 53] and  were successfully fitted to the following expressions:
|B8* (Fitted)||B8* (Simulation)|
Figure 1. Reduced fourth virial coefficient for oblate spheroids versus elongation, . The solid and doted curves are the present work and data obtained by numerical method .
Figure 2. Reduced fourth virial coefficient for prolate spheroids versus elongation, . The solid and doted curves are the present work and data obtained by numerical method .
Figure 3. Reduced fifth virial coefficient for oblate spheroids versus elongation, . The solid and doted curves are the present work and data obtained by numerical method .
Figure 4. Reduced fifth virial coefficient for prolate spheroids versus elongation, . The solid and doted curves are the present work and data obtained by numerical method .
Figure 5. Reduced sixth virial coefficient for oblate spheroids versus elongation, . The solid and doted curves are the present work and data obtained from Monte Carlo integration method .
Figure 6. Reduced sixth virial coefficient for prolate spheroids versus elongation, . The solid and doted curves are the present work and data obtained from Monte Carlo integration method .
Figure 7. Reduced seventh virial coefficient for oblate spheroids versus elongation, . The solid and doted curves are the present work and data obtained from Monte Carlo integration method .
Figure 8. Reduced seventh virial coefficient for prolate spheroids versus elongation, . The solid and doted curves are the present work and data obtained from Monte Carlo integration method .
Table 4. Reduced virial coefficients for hard sphere as limiting case of hard ellipsoid obtained at , in Eqs. (41)- (43) and at in Eqs. (44)- (48). The exact values of hard sphere are obtained from Ref. .
|Reduced virial coefficients||Bn* of hard phere from proposed expressions (41)-(43)||Bn* of hard sphere from proposed expressions (44)-(48)||Exact Bn* for hard sphere|
In Figures. 9-10 the results are compared with the Monte Carlo data. In Table 4 the virial coefficients of hard sphere obtained from Eqs. (41)- (43) as and Eqs. (44)- (48) as , are compared with the exact value of Ref. . As it is being seen, the results are in good agreement. The usefulness of Eqs. (41)- (48) is that they give the permission for estimating the virial coefficients of hard prolate and oblate ellipsoids up to , even for values of and for which no numerical estimate is available. We estimated that the Eqs. (41)- (43) were useful for some convex molecules with .
4. Equation of State of Hard Ellipsoid Isotropic Liquid
4.1. Importance of the Virial Coefficients in the Equation of State of Ellipsoidal Molecular Liquid
More recently, interests in studying the equation of state of non – spherical hard body fluid have increased. The availability of virial coefficients of hard ellipsoids together with the computer simulation data calculated by Frenkel and Mulder  makes it possible to test several types of analytical equation of state for hard convex body fluids. Various equations of state have been proposed in which the compressibility factor is expressed in terms only of the packing fraction, and the shape factor, [19-21]. Such approaches by using a single shape factor must be ultimately inadequate, though for systems of small molecular anisotropy it may be sufficiently accurate. The modified scheme of Naumann and coworkers  might be more successful generally when extreme shapes must be considered. These authors have included a further shape parameter that distinguishes the differences between prolate and oblate molecules.
In this subsection, the importance of virial coefficients are shown in the equation of state of hard ellipsoid fluid and authors’ proposed equations are presented. The virial compressibility factor can be written as follows:
Eq. (49) is an infinite expansion and when higher orders are considered, compressibility factor may get more accurate and closer to the simulation data. We used the virial coefficients up to the eighth as determined in this work. As shown in Figures. 11 to 16, contribution of higher virial coefficients exceeds the accuracy of the compressibility of hard ellipsoid fluids. The obtained compressibility factors at level theory have been compared by the simulation data of Ref. . The convergence of level theory in section V will be discussed.
Figure 11. Compressibility factor of isotropic hard ellipsoid fluid versus the volume fraction of oblate hard ellipsoids with and different virial coefficients contribution. The Monte Carlo simulation data are from Ref. 
Figure 12. As in Figure 11 but for.
Figure 13. As in Figure 11 but for .
Figure 14. As in Figure 11 but for prolate hard ellipsoids with.
Figure 15. As in Figure 11 but for prolate hard ellipsoids with .
Figure 16. As in Figure 11 but for prolate hard ellipsoids with .
4.2. Optimization of the Equation of State of Hard Ellipsoid Fluid
More accurate compressibility factor may be obtained by considering the higher order virial coefficients, as are shown in Figures. 11 to 16. There are more difficulties in calculation of higher order of coefficients. Thus, a limited number of coefficients can be used and the suitable equation of state be deduced. Vega  proposed the following equation:
By adding calculated up to , the following equation have been used:
In addition, we have proposed the equation:
Eq. (51) which is shown in Figures. 11 to 16 is the virial expansion truncated up to eighth coefficient. In Eqs. (51) and (52), we would consider the first eight terms of virial series as determined in this work. In Eq. (52), contribution of virial coefficients higher than the eighth coefficient are obtained by subtracting from the Carnahan-Starling equation the contribution of the first eight virial coefficients of hard spheres. In Figures. 17 to 22, results of the expressions (50) to (52) are plotted and compared with simulation data of the Ref. . The predictions obtained by our new equations of state are more suitable than Eq. (50).
Figure 17. Comparison of the compressibility factor of isotropic oblate hard ellipsoid fluid obtained from introduced equations (50)- (52) and simulation results , for . Eqs. (50) to (52) are due to Vega, our and our proposed expression, respectively.
Figure 18. As in Figure 17 but for .
Figure 19. As in Figure 17 but for .
Figure 20. As in Figure 17 but for prolate hard ellipsoid fluid with
Figure 21. As in Figure 17 but for prolate hard ellipsoid fluid with .
Figure 22. As in Figure 17 but for prolate hard ellipsoid fluid with .
As it is seen in Figure. 23, there are discrepancy with simulation results at high elongation and packing fractions. Therefore, we have proposed another equation of state for hard ellipsoidal molecules. A weighted average of Eqs. (50) and (51) gives:
The indices (50) and (51) are concerned Vega and our equations of state, respectively. We substituted the optimum value of in the proposed expression and are given by:
In Figures. 23 to 25, the new EOS has been compared with the Eqs. (50), (51) and Simulation data in Ref. .
Figure 23. Comparison of the new compressibility factor of isotropic prolate hard ellipsoid fluid with Eqs. (50), (51) and available simulation data  for .
Figure 24. As in Figure. 23 but for.
Figure 25. As in Figure. 23 but for.
5. Results and Discussions
The proposed expressions for sixth to eighth reduced virial coefficients based on and , have shown good agreement with available data obtained by Monte Carlo integration method. This agreement has been indicated in Tables 1- 3. Since the calculation of and need the first and second kinds of elliptic integrals, the sixth- eighth virial coefficients of hard convex molecules can be calculated analytically. In the case of prolate and oblate hard ellipsoids, the viral coefficients of hard ellipsoids versus an elongation parameter, , can be calculated.
The reduced virial coefficients from fourth up to eighth, in terms of elongation are shown in Figures. 1 to 10. The results are in agreement with data of Refs. [31-33, 53] obtained by Monte Carlo method for . Available data to obtain a reasonable expression for eighth virial coefficient were limited. The results of approximate expressions (43) and (48) are shown in Figures. 9-10 and Table 3, for and a good agreement with Monte Carlo data of Refs. [32, 53] can be seen. In Table 4, the virial coefficients for a hard sphere as obtained from limiting case of our proposed equations are compared with the values of Ref. . The results are in good agreement. Also, converged to exact values of hard sphere at , in Eqs. (41)- (43) and at in Eqs. (44)- (48).
We could not find any virial coefficient data for large hard ellipsoids of elongations in literature to compare with our results. Hence, the validity of extrapolation of Eqs. (41)- (48) cannot be guaranteed.
Virial coefficients play an important role in calculation the equation of state of hard convex fluids. In this study, the convergence of the virial series, and the importance of the high-order virial coefficients to test the proposed equations of state were investigated. As it is shown in Figures. 11 to 16, accuracy of EOS of hard prolate and oblate ellipsoidal fluid depends on order of virial coefficients. For aspect ratios , the results converged toward the simulation results at level and reached at and , respectively. For , the results converged toward the simulation results at level; however, approximately reached them at level. The virial expansions up to the eighth order, Eq. (51), are convenient equation of state for hard oblate and small prolate molecules, approximately, as they are shown in Figures. 17 to 22. The higher virials have only a small effect on the equation of state of oblate molecules. However, this is not true in the case of large prolate molecules, and virial expansion up to eighth term is not appropriate equation of state. In order to achieve convergence, we have suggested Eq. (54). Figures. 23 to 25 show the accuracy of Eq. (54) at medium elongations as it is compared with simulated data.
This new EOS provides a better agreement than the other expressions at higher elongations and densities; however, there is failure at . For the extreme aspect ratio , as shown in Figure. 16, the level of viral series gives better agreement than other proposed equations of state with simulations. The new proposed EOS can be applied to hard biaxial ellipsoid and sphero-cylinder isotropic fluids. Also, equation of state of ellipsoid nematic fluid can be calculated.
The complexity of calculations for high order virial coefficients of ellipsoids makes it difficult to obtain accurate analytical high coefficients and equation of state for such systems. In this study, the virial coefficients up to third order were calculated by using analytical method . For higher ones, the numerical values were taken from available data based on Monte Carlo integration method in Refs. [31- 33, 53, 58]. The approximate two and one shape parameter analytical expressions of the hard convex and uniaxial hard ellipsoid molecules have been obtained, respectively, by fitting the available virial coefficients data obtained numerically, up to eighth order. The two set of approximate expressions were examined for prolate and oblate hard spheroids with aspect ratios .
Moreover, the analytical expressions for the equation of state of hard ellipsoid fluids were proposed. For isotropic hard spheroids up to medium elongations, the new proposed equation of state was in good agreement with the simulations. For large elongation, , our new EOS failed to converge and truncated virial series was close to the simulation results. In addition, our proposed equation of state showed a better agreement as compared to others. This new equation was used for both prolate and oblate ellipsoid fluids.
The authors would like to thank the Yasouj University Research Council for partial financial support. This work was supported by the Yasouj University Research Council under grant [Gryu-8913110].