Phase equilibria and interfacial properties of two-dimensional Yukawa fluids

Molecular dynamics simulations of two-dimensional soft Yukawa fluids are performed to analyze the effect that the range of interaction has on coexisting densities and line tension. The attractive one-component fluid and equimolar mixtures containing positive and negative particles are studied at different temperatures to locate the region where the vapor-solid and vapor-liquid phases are stable. When the range of interaction decreases, the critical temperature of the attractive one-component systems decreases. However, for the charged mixtures it increases, and this opposite behaviour is understood in terms of the repulsive interactions which are dominant for these systems. The stable phase diagram of two-dimensional fluids is defined for smaller values of the decay parameter \lambda\ than that of fluids in three dimensions. The two-dimensional attractive one-component fluid has stable liquid-vapor phase diagram for values of \lambda<3, in contrast to the three-dimensional case, where stability has been observed even for values of \lambda<15. The same trend is observed in equimolar mixtures of particles carrying opposite charges.


Introduction
The location of stable phases in systems with molecular interactions has been a subject of great interest in many applications such as simple liquids, electrolyte solutions, colloidal suspension and water. It is well known that fluids with only repulsive interactions have a phase diagram that shows fluid-solid transition [1]. Attractive forces are needed to develop a vapor-liquid phase separation [2]. Several methods, among them, molecular simulations, theory of liquids and experiments are used to obtain information about the phase diagram of atomic and molecular fluids. Although it is possible nowadays to obtain the phase boundaries using atomistic simulations of molecules including internal degrees of freedom [3,4], from theoretical and simulation points of view, it is more convenient to use simple potential models such as square well (SW), Lennard-Jones (LJ), Yukawa and Coulomb interactions. The effect that parameters of the potential has on the phase diagram in general, and in particular, on the critical properties can be systematically analyzed by computer simulations. The SW, LJ and attractive hard-core Yukawa (AHCY) potentials have a short-ranged repulsion and the attraction is long or variable-ranged. For these models it has been established in three-dimensional (3D) systems that the critical temperature decreases as the attraction decreases [5][6][7][8][9]. For very small ranges of attraction, when the vapor-liquid or liquid-solid phase diagrams for the AHCY model become metastable, it has been found that the critical point is below the triple point [8].
The phase diagram and interfacial properties of ionic fluids, where the interaction has attraction and repulsion, has also been studied in 3D for the restricted primitive model (an electroneutral equimolar mixture of particles having the same size and carrying opposite charges) [10,11] and for ionic fluids with asymmetry in size and charge [12][13][14]. The effect of range of interaction on the phase diagram of ionic fluids has been analyzed using the Yukawa potential. Mier-y-Terán et al. [15] used the mean spherical approximation theory to calculate the critical properties of mixtures containing positive and negative particles at different range of interaction. They predicted that as the range of interaction decreases the critical temperature increases, it reaches a maximum and then decreases. That finding is counterintuitive in the sense that less attractive potentials increase the critical temperature. The surface tension for mixtures of charged particles using the Yukawa model have also been reported using molecular dynamics [16]. Fortini et al. [17] performed computer simulations to locate the regions where the vapor-liquid and liquid-solid were stable. Their conclusions were that stable liquid-vapor phase diagrams were found for values of the screening parameter λ < 4, in disagreement with the results of Caballero et al. [18] where the phase diagram was stable for λ < 10. The results of Caballero et al. were obtained for a slightly different potential, but according to Fortini et al., the difference in the results were not related with the potential. Hynninen and Panagiotopoulos [19] found that vapor-liquid phase transition of highly charged colloids is metastable with respect to the vapor-solid phase diagram because at high temperatures the interaction becomes purely repulsive.
Apart from the work in 3D, it is also interesting to understand the phase behavior of ionic particles in two-dimensions. Experimental results on quasi two-dimensional (2D) colloidal suspensions show interesting properties [20,21]. The phase diagrams of 2D systems obtained by computer simulations are scarce. The vapor-liquid and fluid-solid phase diagrams of the LJ fluid in two dimensions were obtained by Barker et al. [22] in 1981 using a liquid-state perturbation theory and Monte Carlo simulations. They found that the phase diagram was qualitatively similar to the 3D system. Santra et al. [23] in 2008 studied the nucleation rate of a liquid phase. They used the LJ model and Monte Carlo simulations to validate the predictions of the classical nucleation theory. They calculated for the first time the line tension of the 2D LJ fluid. Later on, Santra and Bagchi [24] obtained the vapor-liquid phase diagram and line tension of the LJ model at different temperatures. The critical temperature from simulations for the 3D LJ is around 2.5 times greater than its value in 2D. The restricted primitive model of ions in two dimensions was studied by Weis et al. [25] using Monte Carlo simulations and integral equations theory. For this system the critical temperature in the 3D model is around 1.3 times its value in 2D. Analyzing the stability of phase diagrams in 2D might be useful to understand the phase separation in 3D and to investigate new phenomena which are not found in 3D.
The main goal of this work is to analyze the effect that the range of interaction has on coexisting densities and line tension of 2D fluids that interact with the soft Yukawa model. To our best knowledge, there have neither been reported any phase diagrams nor surface tension for these fluids. The one-component systems with attractive interactions and the two-component mixture of positive and negative particles are studied at different ranges of interaction.
This work is organized as follows: The potential model and definition of the calculated properties are presented in section 2. Results are discussed in section 3 and finally Concluding remarks and References are given.

Potential model and calculated properties
The soft Yukawa potential is used to simulate the two-dimensional fluids of pure attractive spheres and equimolar mixtures of particles carrying the opposite charge, where λ −1 is a measure of the range of interaction in dimensionless units and q α is the charge of particle α in the mixture case. The particles are all the same size σ. For the attractive one-component fluid q α q β = −1. The factor f min = 1.075 is included to have a potential which is zero at r = σ and close to −1 at the minimum for the attractive pairs. The short-ranged soft model, instead of the hard sphere model, has been used earlier to study interfacial properties of the restricted primitive model [11,26] and ions with asymmetry in size and charge [27]. The results for the soft model were found to be in good agreement with those where the short range repulsion is given by the hard sphere potential. The advantage of the soft model is that it is straightforward to use molecular dynamics of continuous models.
Reduced units are used in this work for the two-dimensional systems: distance r * = r /σ, energy u * = u/ǫ (where ǫ = U (r min )), temperature T * = kT /ǫ, density ρ * = ρσ 2 , time ∆t * = ∆t (ǫ/mσ 2 ) 1/2 , pressure p * = pσ 2 /ǫ and line tension γ * = γσ/ǫ. The density profile, ρ(x), was obtained as, where 〈N (x, x + ∆x)〉 is the average number of particles with position between x and x + ∆x and ∆A is the area of a slab. The line tension of a planar interface, using the mechanical definition of the atomic pressure, [29] is where L x is the length of the simulation cell in the longest direction and P αα (α = x, y) are diagonal components of the pressure tensor. The factor 0.5 outside the squared brackets takes into account the two symmetrical interfaces in the simulation.
The component P xx of the pressure tensor was calculated as, where v xi and m i are the velocity in the x direction and the mass of particle i , respectively, A is the area of the system and r i j = r i − r j with r i being the position of particle i . A similar expression for P y y was used. The force between particles i and particle j is, (2.5)

Results
Extensive molecular dynamics simulations, with a parallel program, were performed on particles interacting with the soft Yukawa potential at different ranges of interaction λ −1 . All the simulations to study inhomogeneous systems were carried out in non-squared simulation cells, keeping the total density and temperature (NAT) constant.  The density profiles for λ = 4 are shown in figure 2. There was not observed a vapor-liquid phase separation but vapor-solid equilibrium in a very narrow range of temperatures, from 0.30 to 0.315. Large density fluctuations were found at a reduced temperature of 0.32 but a well defined vapor-liquid interface was not stabilized. Homogeneous fluids were found for reduced temperatures above 0.35. The longest size of the simulation cell was increased in order to check if the phase separation from vapor-solid to vapor-liquid was not related to finite size effects. The two-dimensional attractive soft Yukawa model seems to have a stable vapor-liquid equilibrium for values of λ around 3, this is contrary to the 3D case where the same phase equilibrium is stable for values of λ less than 15 [31].  The phase diagram of the attractive soft Yukawa potential is shown in figure 3 for different values of λ. As expected, the critical temperature decreases as the range of interaction decreases. At low temperatures the vapor-solid equilibrium is observed for values of λ = 1.8 and 4. The vapor-liquid transition was not found for λ = 4, it is metastable with respect to the vapor-solid equilibrium. The critical density and temperature for λ = 1 and λ = 1.8 were estimated by fitting the coexisting densities to a rectilinear diameter law with a critical exponent of 1/8 according to a two-dimensional Ising model system [32,33] and experimental results for methane [34]. The final results were (0.36,0.72) and (0.38,0.50), respectively. The critical temperature using the coexisting densities for λ = 1.0 might be around T * = 0.8 but the estimated value using the critical exponent of 1/8 is much smaller. The 2D soft Yukawa with λ = 1.8 and LJ models give nearly the same results as shown in figure 3. The line tension of the attractive soft Yukawa model is shown in figure 4 as a function of temperature, its decay almost following a linear function. The critical temperature can be obtained when line tension is zero, the estimated values for λ = 1.0 and λ = 1.8 are 0.53 and 0.83, respectively. The result for λ = 1 is quite different from the value obtained using the coexisting densities and critical exponent of 1/8. The results for the two-dimensional LJ model, also shown in the figure, are in good agreement with the attractive soft Yukawa model with λ = 1.8.

One-component systems
The phase diagram and line tension of 2D and 3D attractive soft Yukawa fluids are completely different. The vapor-liquid phase diagram in 3D is stable for values of λ < 15, i.e., for very short ranges of interaction, whereas in the 2D case the stability is found for λ < 3. Clearly, the strength of the global attraction needed in 2D to produce the separation is quite different from that in 3D.

Mixture of particles carrying opposite charges
Molecular dynamics simulations on equimolar binary mixtures of equal-size particles carrying opposite charges are carried out to analyze the effect that the range of interaction has on phase stability and line tension. The systems contained 1000 particles in a non-squared simulation cell of L x = 150σ and L y = 38σ dimensions. The potential was truncated at 15σ and the reduced time step was 0.0005. The simulation protocol to obtain the coexisting densities and line tension was the same as that used for the attractive soft Yukawa model described above. The systems were followed for at least 100 blocks of 10 6 time steps. The average properties were obtained from the last 40 blocks. The phase diagram for the same system but in 3D has been calculated [17] and it was found that the vapor-liquid equilibrium was stable for values of λ < 4. In order to find the region where the vapor-liquid transition occurred, several simulations were performed in this work using λ = 2. The difference in critical temperature for the restricted primitive model in 2D and 3D was less than 25 % (see the Introduction). We expected to find the same trend in the binary mixtures using the soft Yukawa model. However, we did not find any evidence of vapor-liquid phase separation using the direct simulation of interfaces. A vaporsolid equilibrium was observed for reduced temperatures as low as 0.106. The amount of vapor increases when the reduced temperature rises to 0.112 but the solid phase still remains stable. At T * = 0.116 the particles behave as a homogenous fluid, see figure 5. As in the 3D mixture of charged particles with opposite sign interacting with the hard core Yukawa model, the critical temperature increases as the range of interaction decreases [15,17], i.e., the system increases its liquid density as λ increases. In fact, for values 23002-5 of λ > 4, the vapor-liquid phase diagram is meta-stable with respect to the vapor-solid equilibrium. In the same way, the possibility of finding a vapor-liquid phase transition in 2D charged mixtures would be in the direction of increasing the range of interaction. Therefore, MD simulations were performed using λ = 1.0. In this case, a vapor-liquid phase separation was found in a very narrow range of reduced temperatures, from 0.078 to 0.08. The liquid contains alternated particles with opposite charges and large voids are observed, see figure 6 for T * = 0.79. In the vapor phase, the particles contain a large cluster and some particles form linear chains and rings.  The line tension for binary mixtures was also calculated and the results are shown in figure 8. As found in simple fluids, the line tension decays with temperature.

Concluding remarks
The main conclusion found in this work is that the stability of vapor-liquid phase diagram in 2D fluids requires longer ranges of interactions than in 3D systems in both the one-component soft Yukawa model and the binary mixture of soft spheres carrying charges with the opposite sign. The range of interaction for these models increases when the inverse screening parameter λ decreases. The metastable vaporliquid transition, with respect to the vapor-solid equilibrium, for the one component model in 2D is found for values of λ > 4 while in the 3D case it is found for λ > 15, i.e., for very short ranges of interaction. For the binary mixture in 2D, metastability occurs for λ > 1 in contrast to λ > 4 observed in the 3D case. A possible explanation might be given using the classical nucleation theory [35] arguments. The formation of a liquid phase in that theory requires that particles in a metastable vapor phase should nucleate up to a critical size. The number of particles in a nucleus of 2D fluids for the same range of interaction is smaller 23002-6 Two-dimensional Yukawa fluids than in 3D because in 3D the particles are in a sphere while in 2D they are in a circle. The particles in 2D have to interact longer distances than in 3D for the nucleus to reach the critical size.
On the other hand, the vapor-liquid phase diagram and line tension for the one-component attractive soft Yukawa model are in good agreement with those obtained using the LJ in 2D.