SciELO - Scientific Electronic Library Online

 
vol.65 issue2Thermoactivated cavitation induced in water by low power, continuous wave thulium-doped fiber laser author indexsubject indexsearch form
Home Pagealphabetic serial listing  

Services on Demand

Journal

Article

Indicators

Related links

  • Have no similar articlesSimilars in SciELO

Share


Revista mexicana de física

Print version ISSN 0035-001X

Rev. mex. fis. vol.65 n.2 México Mar./Apr. 2019  Epub Apr 17, 2020

https://doi.org/10.31349/revmexfis.65.190 

Investigación

Time dependent self-diffusion coefficient of methane molecules confined into micropores structure of a sandstone

F. de J. Guevara-Rodrı́gueza 

a Gerencia de Desarrollo de Materiales y Productos Químicos, Instituto Mexicano del Petróleo, Eje Central Lázaro Cárdenas 152, Ciudad de México, 07730, México. e-mail: fguevara@imp.mx


Abstract

In this work, the effect of pore structure of a sandstone on the molecular displacement of confined methane gas is analyzed. It was found that the self-diffusion coefficient of a methane molecule depends on the pore size distribution. In particular, the time dependent self-diffusion coefficient exhibits a maximum which is correlated with the effect of the molecular confinement. It was found that a sandstone with small pores (whose diameter is less than 20Å) traps the gas more efficiently than other sandstones.

Keywords: Porous materials; time dependent self-diffusion; confinement effect

PACS: 47.56.+r

1. Introduction

The pore systems in a sandstone could be one of four basic porosity types 1, namely: (1) intergranular, (2) porosity, (3) dissolution, and (4) fracture. In the first case, intergranular exists as space between detrital grains. In the second case, porosity exists as small pores (less than 2 μm) commonly associated with detrital and authigenic clay minerals. In the third case, dissolution is the pore space formed from the partial to complete dissolution of framework grains and/or cements. Finally, in the fourth case, fracture is the void space associated with natural fractures.

Natural gas is a gaseous mixture containing at least 75 vol.% of methane 2, and initially occupies between 30-90 % of total space of pores in the sandstone in a fresh reservoir. As an example, natural gas initially occupies between 45-60 % of pores volume at the Burgos province at the north of Mexico. Part of the volume is initially occupied by the gas and the rest of the volume is occupied by water at the bottom of reservoir. However, the reservoir is invaded by aquifer water with the gas production on the well. At the end of the productive live of the gas well, Residual Trapped Gas Saturation (RTGS) measure is a key factor to evaluate the additonal gas recovery from a drained gas reservoir, however, the measurements of RTGS exhibit values which are scattered from 5 % to 85 %. This fact represents one of the main uncertainties in the recoverable reserves of the field. To understand the measurements of RTGS some hypotheses are laid out to explain this phenomenon 3, 4. One of them is that during the gas production, water invades into the gas-saturated zone trapping a certain amount of gas independently of pores structure of the sandstone 3. Another hypothesis is that RTGS must decrease for sandstones with high porosity but it is not clear at all because for another sandstones (with similar porosity as the first one), the RTGS increases 4. At this point, the last comment suggests that the pores structure of the sandstone is the key to understand the measure of the RTGS, because at molecular level, the pores size distribution affects the displacement of methane molecules trapped into the intergranular pores or micropores 5, 6. Thus, dynamic properties (for example, the time dependent self-diffusion coefficient of methane molecules) must tell us something about the molecular confinement. In the literature, few manuscripts are focused on the self-diffusion coefficient 7-10 where the diameter of the pores is in the range of 100-500 Å. In this work, the self-diffusion coefficient is analyzed in the case of a confined gas in intergranular pores where their diameter is less or around of 20 Å. This is the goal of the present work.

The paper is divided into three main sections. Section 1 is focused on the construction of the model of the porous material and its characterization through the intergranular pores size distribution. Section 2 is focused on the analysis of the self-diffusion coefficient of the methane molecule confined into the intergranular pores of the porous material of previous section. Finally, conclusions are in the last section.

2. Model of a porous material

In this work, we focused only on intergranular porosity of a material. A micropore size definition, as the pore whose diameter is less than 20Å, was emitted by the commission on colloid and surface chemistry including catalysis of the International Union of Pure and Applied Chemistry (IUPAC) 11. On the other hand, the diameter of a methane molecule is around of σm3.73 Å 12-14 and this value is 18.65 % of the diameter of a micropore in the IUPAC’s definition. In this work, the diameter σ0=18.65 Å of a hypothetical micropore (nears to the value in the IUPAC’s definition) is only used as a unit length along of the rest of the manuscript. Thus, a unit length is 5× the diameter of a methane molecule.

A porous material is modeled as a mixture of hard spheres. In particular, the pore size distribution is analyzed. All models have the same value for the volume fraction occupied by spheres, namely, η=0.6. In this point, the pores structure is related to the number and the concentration of species of hard spheres in spite of the same volume fraction occupied by spheres. In this work, the models emulate a sandstone with high porosity [15,16].

Four models are constructed from a binary mixture of hard spheres and their configuration of species are in Table I. In the same way, another four models are constructed from a ternary mixture of hard spheres, and their specifications are in Table II. The pore size distribution is a function of the configuration of species in the mixture of hard spheres, but we can not speak about it without a previous definition of a pore, i.e. what is a pore? Moreover, if a reasonable definition of a pore is established, then what is its volume? The answer for both questions are in the following section.

TABLE I Configuration of a binary mixture of hard spheres for the model of a porous material. The value of σ0 = 18:65 °A is used as a unit length. 

Specie 1 Specie 2
Model σ1/ σ0 N 1 σ2/ σ0 N 2
A 1.000 500 0.100 1000
B 0.780 500 0.320 1000
C 0.550 500 0.550 1000
D 0.278 500 0.822 1000

TABLE II Configuration of a ternary mixture of hard spheres for the model of a porous material. The value of σ0 = 18:65 Å is used as a unit length. 

Specie 1 Specie 2 Specie 3
Model σ1/ σ0 N 1 σ2/ σ0 N 2 σ3/ σ0 N 3
E 1.000 100 0.500 900 0.250 0
F 1.000 100 0.500 900 0.250 500
G 1.000 100 0.500 900 0.250 1500
H 1.000 100 0.500 900 0.250 2500

2.1. Volume of a pore

The definition of a pore is illustrated in Fig. 1. In this case, the tetrahedron is formed by the centers of four neighbor spheres. Thus, the shape of the pore corresponds to tetrahedron volume without the partial volume of each sphere. The algorithm to calculate the pore size is now described:

  • 1. A sphere in the matrix of hard spheres is selected (and is labeled as the sphere 1). Other three spheres, which are more close to sphere 1, are selected too. The centers of the four spheres are the corners of the tetrahedron (see the Fig. 1), and the tetrahedron volume is calculated with equation

V0=16|r21(r31×r41)|, (1)

  • where the vectors r21, r31, and r41 are defined in the Fig. 1.

Figure 1 Tetrahedron formed by 3 of the more close spheres to the sphere 1

  • 2. In this step a system of coordinate axes is chosen so that the center of the sphere 1 is at the origin, meanwhile the unit vector r^21 is on the axis z, the unit vector r^31 is on the plane xz, and the unit vector r^41 together with r^21 define a second plane that forms the angle Φ with the plane xz. The system of coordinate axes is illustrated in Fig. 2 where the angles θ1, θ2 and ϕ are determined from the following vector operations

Figure 2 Partial volume of a sphere. The center of the sphere 1 is located at the origin. 

r^21r^31=cos(θ1); (2a)

r^21r^41=cos(θ2); (2b)

r^21(r^31×r^41)=sin(θ1)sin(θ2)sin(ϕ). (2c)

The size of the partial volume of sphere 1 is calculated with the formula

ΔV1σ1324[(1-cos(θ2))ϕ+(cos(θ2)-cos(θ1))×arctan(sin(θ2)sin(ϕ)sin(θ2)cos(ϕ)+sin(θ1))], (3)

  • where σ1 is diameter of the sphere 1. For some geometries, equation () provides us the exact formula for the partial volume of a sphere, but in other cases the expression can be used to calculate the value of ΔV1 in an approximately way. Of course, equation () can be also used to compute the partial volume of the other spheres in the tetrahedron by considering, for example, that sphere 2 now plays de role of sphere 1 and so on for spheres 3 and 4.

  • 3. In this final step the size of the pore is calculated with V=V0-iΔVi.

This algorithm is applied on each sphere in the mixture of hard spheres (sandstone model) and the set of values of V will be used to construct the pore size distribution as it is described in the next section.

2.2. Pore size distribution

In order to explore different regions of the sandstone model, a set composed by N=10000 configurations is constructed for any of the models in Tables I and II by using the Monte Carlo algorithm 17. The volume V of the pore is calculated by using the procedure of Sec. 2 for each site in the matrix of the model and for all configurations. The nth element of the histogram (hn) is increased by 1 if the pore volume fulfills with Vn-1<VVn, where Vn=n×10-5Vbox for n=1,2,,NV with NV=1000 and Vbox is the volume of the simulation box. This is the method to construct the histogram with the pore size distribution.

Once the construction of histogram h(V)is ended, in the next step, the histogram is normalized by using the next formula

f(σ)=h(V)A, (4)

where σ is the diameter of a sphere with the same volume of the pore (V), thus the “pore diameter” is related to the pore volume through V=πσ3/6, meanwhile, the denominator (A) is the normalization factor which is defined by

A=0h(V)dσ=29π30h(V)V2/3dV. (5)

In Fig. 3 the pore size distribution as a function of the pore diameter is plotted for the models A, B, C, and D reported in Table I. Model C is composed by a single specie of spheres and, in Fig. 3, its curve shows a single maximum with a narrow distribution around it. This fact tells us that the model C is characterized by a pore set of similar sizes. On the other hand, pore size distribution of models A and B show two maximums which are related to two different main sets of pores but in a narrow distribution as model C. Furthermore, the pore size distributions of models A, B, and C have a similar range of the pore size, namely, σ/σ0(0,1). In contrast, the pore size distribution of model D exhibits a maximum in a range where a pore diameter is bigger than the pore diameter in the other models A, B, and C.

Figure 3.Pore size distribution as a function of the diameter σ of a sphere with the same pore volume. The curves correspond to the models A, B, C, and D reported in Table I.  

The pore size distribution for the models E, F, G, and H (see Table II) are in Fig. 4. In these models, the number of small spheres of the third class increase from model E to model H, meanwhile the configuration of other two species does not change. Thus, the pore size decreases with the number of spheres of the third class and this fact can be to see in Fig. 4 with the shift to the left of the pore size distribution. On the other hand, model E is a binary mixture of hard spheres because the third specie is not present. In this case, the pore size distribution exhibits two maximums related to two main set of pores in the system. Models F, G, and Hare ternary mixtures of hard spheres and, in particular, the pore size distribution of models G and H is like a bell but with a soft ripple. This fact tells us that these pore size distributions have a more complex structure where, for example, model H clearly exhibits three peaks which are related to three main sets of small pores (if these pores are compared with the pore size in other models, namely, E, F, and G).

Figure 4 Pore size distribution as a function of the diameter σ of a sphere with the same pore volume. The curves correspond to the models E, F, G, and H reported in Table II

Another relevant effort to characterize the porous material is by using a tracer for which its molecular displacement is affected by the temperature, the density of the fluid, and the properties of the sandstone, i.e., the pores size distribution (among other properties like the porous connectivity which is not analyzed in this work). In the next section, the model of a methane gas confined into the micropores of the sandstone is analyzed and, in particular, a dynamic property of the methane molecules will be used to characterize the system, namely, the time dependent self-diffusion coefficient of methane.

3. Sandstone and methane gas model

The construction of a sandstone model was discussed in the previous Sec. 2 and the pore size distribution was used to characterize it. At this point, the hard sphere condition for the sites of the sandstone model was only used to construct the porous material. In this section, and for the rest of the work, the hard sphere condition is removed and is substituted by a Lennard-Jones site (where its diameter is the same of previuos hard sphere). From a static matrix of sites the methane molecules (modeled with spheres) are initially placed in it in a random way. In the next step, all overlaping spheres are moved until all overlaps are completely removed. In this way, the model of a sandstone and confined gas is initially constructed. In this section, the time dependent self-diffusion coefficient of methane molecule of the gas, which is confined into the micropores of the sandstone, is discussed. The full model is illustrated in Fig. 5 where Nm and Ns are the number of methane molecules and the number of sites of the rock model (where, in the way, the sites are static over the time), respectively. All spheres in Fig. 5 are in a cubic simulation box. The sandstone models are reported in Tables I and II, and all of them shared the same value η=0.6 for the fraction of occupied volume, and therefore, the free volume which is available for methane gas is Vf=(1-η)Vbox where Vbox is the total volume of the system and corresponds to the size of the simulation box. On the other hand, as was mentioned in Sec. 1, the range of values of the RTGS measure is from 5 % to 85 % of the free volume (Vf) and in this analysis x=0.05 (5 %) is considered for the RTGS value. Thus, the configuration of all sandstone and methane gas models are in Table III where σm/σ0=0.2 is the methane diameter.

Figure 5 Illustration of the simulation box with the sandstone model (gray spheres) and the methane gas (blue spheres).  

TABLE III Configuration of the sandstone and methane gas models. Ns is the number of sites in the sandstone (where, in the way, the sites are static over the time); Nm is the number of methane molecules; Vbox is the volume of the cubic simulation box; Vs = ηVbox is the volume occupied by the sites of the sandstone; Vm = xVf is the volume occupied by the methane molecules; and Vf = Vbox - Vs is the free volume available for the gas. 

Model Ns Nm Vbox/σ03 Vs/σ03 Vm/σ03 Vf/σ03
A 1500 2088 437.2 262.3 8.7 174.9
B 1500 1125 235.7 141.4 4.7 94.3
C 1500 1040 217.8 130.7 4.4 87.1
D 1500 2359 494.1 296.4 9.9 197.6
E 1000 885 185.4 111.3 3.7 74.2
F 1500 918 192.3 115.4 3.8 76.9
G 2500 983 205.9 123.5 4.1 82.4
H 3500 1048 219.5 131.7 4.4 87.8

The gas molecules interact between them and with the sites of the sandstone. In particular, the total potential energy (U) is calculated and approached with the sum of the pair interactions, i.e.,

U(r1,,rNm)=i=1Nm[12j=1jiNmφ(|ri-rj|)+j=1Nsψ(|ri-Rj|)], (6)

where ri and Rj are the center of the ith methane molecule and the jth sphere in the sandstone, respectively. Moreover, the energy between a pair of molecules is approximately calculated with the Lennard-Jones potential

φ(r)=4εm[(σmr)12-(σmr)6], (7)

where the parameters are εm/kB147.9K and σm3.73Å 18, 19. In the same way, the pair interaction between a gas molecule and a site of the rock is also calculated with the Lennard-Jones potential

ψ(r)=4εR[(σRr)12-(σRr)6], (8)

where the parameters are εR=0.5εm and σR=(σs+σm)/2. For the complete system (see Fig. 5) periodic boundaries are considered and cut-off radius is used to compute the Lennard-Jones potential, i.e., φ(r)=0 if r>2.5σm and ψ(r)=0 if r>2.5σR. Finally, the total force in the ith methane molecule is computed with fi=-iU.

3.1. Molecular dynamics algorithm

Molecular dynamics simulation of confined gas into the porous material is performed by using the reversible in time algorithm which was proposed by Martyna 20-23. In this case, procedure to do the numerical integration of movement equations is described with the following set of equations

vi*=vit0+Δt2mfit0; (9a)

ri(t1)=ri(t0)+Δtvi*; (9b)

vi(t1)=vi*+Δt2mfi(t1), (9c)

where ri(t0), vi(t0), and fi(t0) are the position, the velocity, and the force on the gas particle at time t0. m is the mass of the particle i (where i=1,2,,Nm) and is proportional to the molecular weight of methane: 16.0426 g/mol. At the time t1=t0+Δt, Eqs. () and () enable us to update the position of all molecules from ri(t0) to ri(t1). Once the new configuration is established then the force fi(t1) is computed from the force field mentioned in Sec. 3.1 and the resulting force is used to update the velocity of the gas particles with Eq. (). In this way, movement equations are integrated over the time and Δt=0.001ps is the time step which is used in this work. Furthermore, the Nosé-Hoover thermostat is coupled to the above numerical integrator to preserve the thermal equilibrium 24-25.

After the system is constructed with an initial valid configuration, then the molecular dynamics simulation is performed by generating at least 2000 steps over time starting from its initial configuration to those configurations which correspond to a system in thermal equilibrium at the temperature of T=433.15K (which is a typical temperature of a gas reservoir). Once the system reaches the thermal equilibrium then, in the next step, a dynamic property of the confined gas is calculated trough a second process of its molecular dynamics simulation. The details of the method are in the next section.

3.2. Self-diffusion coefficient

In this step a new molecular dynamics simulation is performed by generating a sequence of 10000 consecutive steps by using the algorithm discussed in previous Sec.3. However, the information of the dynamic state of the system is saved in a external file every 5 steps over time, i.e., the output file contains Nc=2000 configurations of the system.

After the output file has been constructed, in the next step, the time dependent self-diffusion coefficient is calculated (from configurations which are in the output file) by using its formal definition, namely,

D(t)=|r(t)-r(0)|26t, (10)

where D(t) is the mean squared displacement divided by the time and it is a measure of the deviation of the position of a particle with respect to a reference position over time 26-28. Here, is the average in the NVT ensemble. Equation () is equivalent to

D(tn)=1Nmi=1Nm1Nc-n×j=1Nc-n|ri(tj+tn)-ri(tj)|26tn, (11)

where tn=5nΔt and n=1,2,,Nc-1. In this way, time dependent self-diffusion coefficient can be used as a measure of the portion of the porous material which is “explored” by the gas particles.

Self-diffusion coefficient of confined gas into the pores of a binary mixture of spheres (sandstone model) are in Fig. 6.

Figure 6 Self-diffusion coefficient of confined methane molecules into the porous material. The sandstone model corresponds to a binary mixture of spheres which its configuration is in Table I.  

Curves plotted in Fig. 6 correspond to the cases A, B, C, and D which are reported in Tables I and III. Clearly, self-diffusion coefficient depends on the structure of pores as it is seen in Fig. 6. Furthermore, the sandstone models A and D have the mayor size of free volume which is available to the gas particles (see Table III), thus, self-diffusion coefficient increases apparently with the free volume. The exception is the case D. The size distribution of pores in model D (see Fig. 3) exhibits a main set of pores, meanwhile the size distribution of pores in model A exhibits two set of pores. In model A, the pores are small with respect to the pores in model D but the two set of pores in model A have a clear overlapping where, perhaps, its effect is that the molecular displacement has the best performance with respect to the other models. The most relevant feature of the self-diffusion coefficient curves, which are plotted in Fig. 6, is a notorious peak at intermediate times. If a particle trapped in a pore cannot escape from it immediately then its displacement is around of the geometrical center of the pore and the time in it increases until the particle finally scape. With this picture in mind, a narrow peak, but a high peak with respect to the long-time value of the self-diffusion coefficient, is correlated with a methane molecule confined into the pore over the time. Thus, the self-diffusion coefficient of model C corresponds to a typical case of a gas trapped in an efficiently way by the sandstone. In contrast, model A corresponds to a sandstone where the two set of main pores facilitate the methane molecules displacement in the pores structure. Other feature is the position of the peak on the time. In this case, if the pores are small in the pores structure then the position of the peak is found quickly at the early-times, i.e., the gas particles “detect” the confinement more quickly if they are into small pores.

The scenario in the other sandstones, which are modeled with a ternary mixture of spheres, is also similar to the binary cases and the self-diffusion coefficient curves of the methane molecule, are plotted in Fig. 7. All the features discussed previously for the binary models of a sandstone are present in the ternary models. However, the curves in the ternary models are in general below with respect to the binary mixture cases, thus the size of pore in the ternary mixture of spheres are small with respect of the a binary mixture of spheres. Furthermore, the self-diffusion coefficient of models E, F, G, and H exhibits a high narrow peak in all curves which are plotted in Fig. 7, thus the sandstone (which is modeled with a ternary mixture of spheres) traps the gas more efficiently than the binary mixture of spheres.

Figure 7 Self-diffusion coefficient of confined methane molecules into the porous material. The sandstone model corresponds to a ternary mixture of spheres which its configuration is in Table II.  

On the other hand, the initial slop of self-diffusion coefficient in all cases reported in Fig. 6 and also in Fig. 7 have the same value. To clarify this point, short-time regime 8 is defined by ξD0t/σ021 where D0 is the molecular self-diffusion coefficient at short-time regime and bulk condition, i.e., D 0 corresponds to molecular self-diffusion of an un-confined gas. Thus, ξ1 means that molecular displacement is extremely less than the pore size and, in this regime and from Eq. (), self-diffusion coefficient is approached with

D(t)kBT2mt, (12)

where kB is the Boltzmann constant, T is the temperature of reservoir, and m is the mass of methane molecule. Clearly, initial slope is a function of temperature which is constant in this analysis.

4. Conclusions

The concentration of methane molecules is Nm/Vf=6x/πσm3 and, therefore, it is constant in all models in Table III. Moreover, temperature is another constant in the analysis of self-diffusion coefficient in Sec. 3.2. Thus, the resulting time dependent self-diffusion coefficients, which are plotted in Figs. 6 and 7, are functions of the pore structure of sandstone. In particular, the sandstones with small pores traps the gas more efficiently than other pores structure. On the other hand, self-diffusion curves in all cases exhibit a maximum and the shape of the peak is correlated to the effect of the molecular confinement into a pore. In fact, if peak is high and narrow then it means that the particle can not to scape from the pore immediately. Thus, the molecular displacement is around of the geometrical center of the pore and the time in it increases until the particle finally scape. At the long-time regime (ξ1) the self-diffusion coefficient approaches to a plateau defined by D(t)DL if t. In this case, the value of D L increases with the size of free volume available to gas. Finally, in this work, the value of RTGS is 5% and it is found that the self-diffusion coefficient depends on the pores structure. However, the inverse sentence could be not valid, i.e., for two sandstones where the self-diffusion coefficients are similar then the RTGS measure depends on the pores structure. This question is still in the air and requires to be analyzed in a future work.

Acknowledgments

The author thanks to Instituto Mexicano del Petróleo for its support through project D.61072.

References

1. Edward D. Pittman, Porosity, Diagenesis and Productive Capability of Sandstone Reservoirs, Aspects of Diagenesis SEPM Society for Sedimentary Geology, (1979). 9781565761568. [ Links ]

2. Marco Tagliabue et al., Chem. Eng. J. 155 (2009) 553-566. [ Links ]

3. Nicola Bona et al., SPE (2014) October doi:10.2118/170765-MS [ Links ]

4. G. Hamon and K. Suzanne and J. Billiotte and V. Trocme, SPE January (2001) doi:10.2118/71524-MS. [ Links ]

5. B. Hosticka and P. M. Norris and J. S. Brenizer and C. E. Daitch, (1998) J. Non-cryst. Solids 225 293-297. [ Links ]

6. Sen Wang Farzam Javadpour and Qihong Feng, Fuel 171 (2016) 74-86. [ Links ]

7. W. D. Dozier, J. M. Drake and J. Klafter, Phys. Rev. Lett. 56 (1986) 197-200. [ Links ]

8. Rustem Valiullin and Vladimir Skirda, J. Chem. Phys. 114 (2001) 452-458. [ Links ]

9. Kourosh Malek and Marc-Olivier Coppens, Phys. Rev. Lett. 87 (2001) 125505. [ Links ]

10. David J. Bergman, Keh-Jim Dunn, Lawrence M. Schwartz and Partha P. Mitra, Phys. Rev. E 51 (1995) 3393-3400. [ Links ]

11. J. Rouquerol et al., Pure & Appl. Chem. 66 (1994) 1739-1758. [ Links ]

12. Mohsen Abbaspour and Hamed Akbarzadeh and Sirous Salemi and Mousarreza Abroodi, Phys. A 462 (2016) 1075-1090. [ Links ]

13. William L. Jorgensen and Jeffry D. Madura and Carol J. Swenson, J. Am. Chem. Soc. 106 (1984) 6638-6646. [ Links ]

14. Marcus G. Martin and J. Ilja Siepmann, J. Phys. Chem. B 102 (1998) 2569-2577. [ Links ]

15. G. E. Archie, Bull. Amer. Assoc. Petrol. Geol. 36 (1952) 278-298. [ Links ]

16. H. J. Sander Bull. Amer. Assoc. Petrol. Geol. 51 (1967) 325-336. [ Links ]

17. Pascal H. Fries and Jean-Pierre Hansen, Mol. Phys. 48 (1983) 891-901. [ Links ]

18. Michael P. Allen and Dominic J. Tildesley, Computer Simulation of Liquids: Second Edition (Oxford University Press 2017). 9780198803195 [ Links ]

19. Daan Frenkel and Berend Smit, Understanding Molecular Simulation: Second Edition (Academic Press 2001). 9780122673511 [ Links ]

20. M. Tuckermar and B. J. Berne and G. J. Martyna, J. Chem. Phys. 97 (1992) 1990-2001 [ Links ]

21. D. Levesque and L. Verlet J. Stat. Phys. 72 (1993) 519. [ Links ]

22. Trond Ingebrigtsen and Ole J. Heilmann and Soren Toxvaerd and Jeppe C. Dyre, J. Chem. Phys. 132 (2010) 154106. [ Links ]

23. William Graham Hoover and Carol Griswold Hoover World Scientific Publishing (2012). 9789814383165 [ Links ]

24. S. Nosé Prog. Theor. Phys. Suppl. 103 (1991) 1-46. [ Links ]

25. W. G. Hoover Phys. Rev. A 31 (1985) 1695-1697 [ Links ]

26. Jean-Pierre Hansen and I.R. McDonald Theory of Simple Liquids Academic Press 9780123870322 (2013). [ Links ]

27. Rajneesh K. Sharma and K. N. Pathak and K. Tankeshwar and S. Ranganathan, Phys. Chem. of Liq. 29 (1995) 59-68. [ Links ]

28. David M. Heyes and Jack G. Powles and Gerald Rickayzen, Mol. Phys. 100 (2002) 595-610. [ Links ]

Received: August 02, 2018; Accepted: October 06, 2018

Creative Commons License This is an open-access article distributed under the terms of the Creative Commons Attribution License