Return

The microcluster simulations using MD method [2] generally require disposition of target microcluster in a sufficiently large size cell. And the cell is applied to free boundary condition. However, in cases where the size of such target microcluster is small and molecules in the microcluster are monatomic, the energy source is limited only to intermolecular interaction and kinetic energy. Under such condition, low temperature preset for the simulations may prevent the microcluster from obtaining sufficient energy and such microcluster may remain unchanged. To avoid the result that such microcluster is not able to obtain sufficient energy from outside, the model [3] disposing numerous career gases in addition to molecules which constitute the microcluster in the cell is often used. And the cell is applied to periodic boundary condition. By this method, when career gas particles collide with microcluster, the microcluster can obtain energy from the exterior. As a result, the phase transition phenomenon is promoted. However, this method may cause the behavior of microcluster with preset pressure and the career gas number. And accurate thermodynamic quantities of microcluster may not be calculated by existence of career gas. Accordingly, since the phase transition phenomenon is intentionally promoted using career gas, the calculation result is not necessarily right.

Therefore in this paper, argon microcluster constituted monatomic was taken up for an example, MD simulations using career gas were calculated, and the thermodynamic quantities of pure Ar microcluster on phase transition were obtained. In regard to career gas, we used helium molecule belonging to rare gases. However, Lennard-Jones parameters of career gas were Ar ones for the purpose of simplifying the model. Specifically, interaction terms of internal energy, enthalpy, heat capacity, entropy and free energy were obtained. We also performed MC simulations [4] without career gas and obtained the quantities like the MD method, in order to compare the results from two methods and evaluate the suitability of MD method.

With respect to temperature, we preset it at 1 K (= 0.00834 e/

In eq. (2),

Table 1. Setting conditions required for MD method.

Numerical Integration | Gear Method |

Temperature Control | Velocity Scaling |

Ensemble | NTP |

Initial Density | 0.00142g/cm^{3} |

MD Steps | 100000 |

Time Step | 1 fs |

3D Periodic Boundary Condition | |

Cubic Cell |

In addition, the software and environment of the personal computer used for MD simulations are shown in the following Table 2.

Table 2. PC spec. for MD simulations.

Software | WinMASPHYC-Pro1.0 |

OS | Windows98SE |

CPU | PentiumIII600MHz |

MEMORY | SDRAM256MB |

The thermodynamic quantities per one Ar molecule were calculated in the following procedure. First, the thermodynamic quantities of Ar and He mixed system were calculated. Next, the thermodynamic quantities of the pure He system were calculated. And we subtracted the latter from the former. Furthermore, the result was divided by number of Ar molecules (see eq.3):

The reliability of the thermodynamic quantities per one Ar molecule solved by using the above procedure is explained in the following section. The thermodynamic quantities were shown separately by number of career gas in Figures 1 - 3.

Figure 1. *U*_{e} values per one Ar molecule varying with number of career gas, as a function of temperature. (*N*=8)

Figure 2. *P*_{e}*V* values per one Ar molecule varying with number of career gas, as a function of temperature. (*N*=8)

Figure 3. *V* values per one Ar molecule varying with number of career gas, as a function of temperature. (*N*=8)

Hereinafter how to calculate thermodynamic quantities of each figure are explained briefly. Concerning the values *U*_{e} and *P*_{e}*V* already shown in eq. (1), the ideal gas energy values were subtracted from the total energy values, because calculation of the interaction terms generated by Ar is the purpose of this research.

For comparison purposes, Figures 1 - 3 include the results obtained by NTV ensemble. In the NTV ensemble, the career gas and the periodic boundary condition are disabled. The velocity scaling method is used for a temperature control. On the result of NTV ensemble, *U*_{e} value increases around *T* = 0.48 and 0.72 e/*k*_{B} when the structure of Ar microcluster changed. However decomposition of Ar microcluster could not be observed within the range of set-up temperature on NTV ensemble.

Next, in NTP ensemble, periodic boundary condition and disposition of career gas are enabled. In this setting, *U*_{e} values per one Ar molecule increase drastically from the point at which the temperature exceeds *T* = 0.3 e/*k*_{B} and such values finally reach *U*_{e} = 0 e. This result means that the microcluster has collapsed and evaporated with first order phase transition. However the number of career gas influences *U*_{e} value. The *U*_{e} curves about the different numbers of career gas show relevance between increase of the number of career gas and drop in the phase transition temperature. Although such number increase does not unlimitedly bring the phase transition temperature down it seems that such temperature converges into the real phase transition temperature. Such tendency is explained by the following observations. In the case of 50 He molecules, the phase transition temperature was around *T* = 0.43 e/*k*_{B} due to obvious lack of career gas number. However in the cases of 100 and 200 He molecules, their respective phase transition temperatures reached *T* = 0.34 and 0.32 e/*k*_{B}, and both temperatures showed little difference. At the same time, in the case of 200 He molecules their phase transition temperatures ranged shortly from *T* = 0.31 to 0.32 e/*k*_{B} whereas on the contrary in case of 100 He molecules the temperatures ranged widely from *T* = 0.34 to 0.45 e/*k*_{B}. It is considered that lack of energy is the cause. Ar microcluster became unstable due to high energy generated, so it released one molecule and became stable again. At the same time the energy that microcluster received from career gas 100 He molecules became insufficient to maintain its phase transition, and such transition terminated at such temperature. From these results, the large number of career gas brings the phase transition temperature down to its accurate temperature, and shortens the phase transition temperature range and the required time. Figure 2 shows almost constant *P*_{e}*V* value notwithstanding change of the number of career gas molecules. According to Figure 3, the assembly volume grows in proportion to the number of career gas. However, the pure Ar system volume was not influenced by the number of career gas molecules.

Figure 4. *U*_{e} values per one Ar molecule varying with pressure intensity, as a function of temperature. (*N*=8)

Figure 5. *P*_{e}*V*values per one Ar molecule varying with pressure intensity, as a function of temperature. (*N*=8)

Figure 6. *V*values per one Ar molecule varying with pressure intensity, as a function of temperature. (*N*=8)

Figures 4 - 6 show the thermodynamic quantities separated by change of pressure. In the figures, we used the system that contains 200 He molecules because the most accurate thermodynamic quantities were obtained by using 200 He molecules as career gas from Figures 1 - 3.

Each thermodynamic quantity was calculated under pressures of 0.1, 1.0 and 10 MPa. According to Figure 4, at 0.1 MPa level, high energy state has been generated from the beginning and movement of Ar molecules is perceived at around *T* = 0.1 e/*k*_{B}. In addition, such energy has been consistently unstable due to low pressure and high amplitude of expansion and contraction of the cell. And the starting point of phase transition has been unascertained because of gradual increase of *U*_{e} value until reaching zero. On the contrary, at the 10 MPa level, high pressure prevented the cell from expanding sufficiently and both Ar microcluster and career gas were forced to fit themselves in the small cell despite the temperature rise. On the preceding condition, even if such Ar had been liberated from Ar microcluster, it collided with career gas and reattached itself to microcluster immediately after such liberation. As a result, the phase transition temperature was considerably high. In accordance with Figures 5, 6, 0.1 and 10 MPa curves are unstable at *P*_{e}*V* value and particularly such instability is outstanding after the phase transition compared to the stability of the pressure 1 MPa curve.

According to Figures 1 - 6, it is ascertained that the 8 Ar molecules microcluster selected for this simulation requires the career gas number to be 200 molecules, pressure to be 1 MPa and initial density of cell to be one hundredth of liquid density 0.0142 g/cm^{3}.

Table 3. Setting conditions of MD method and MC method.

Calculation Technique | MD Method | MC Method |
---|---|---|

Numerical "Integration" | Gear Method | Metropolis Method |

Ensemble | NTP | NTV |

Used Molecule | Ar(8)+He(200) | Ar(8) |

Initial Volume | 3318.03 s^{3} | 3795.36 s^{3} |

Pressure | 1MPa | -- |

Periodic Boundary Condition | Enabled | Enabled |

Time Step | 1fs | -- |

Table 4. PC spec. for MC simulations.

Development language | Compaq Visual Fortran6.5 |

OS | Windows98SE |

CPU | PentiumIII600MHz |

MEMORY | SDRAM256MB |

The development language and spec. of the personal computer used for MC simulations are shown in the Table 4. The personal computer is the same as that used for calculation of the MD method with the section 2. The program for the MC method was developed by us.

In the case of MC method, career gas was not disposed on the contrary to the MD method. Sufficient number of steps in MC simulations leads to accurate results, so it is expected that this calculation gives us thermodynamic quantities of pure Ar. In addition, we used the volume measured immediately after completion of phase transition for using MD method in MC simulations with molecular extension because the volume is constant according to NTV ensemble. Each thermodynamic quantity calculated under such setting conditions using MC method compared to MD method is shown in Figures 7 - 14 below.

In Figure 7, average internal energy *U* obtained by using MD method in two systems; Ar and career gas mixture system and pure career gas system is shown. In Figure 8, the product of pressure and volume *PV* is shown as Figure 7. On one hand, the mixed system curves change with phase transition, on the other hand, the pure career gas system curves are linear. As referring to eq. (3), we obtain thermodynamic quantity per one Ar molecule by subtracting pure career gas system value from mixed system value described in these figures and dividing it by Ar molecules number 8. It is not easy to determine whether the bold procedure using subtraction gives an accurate calculation result of Ar value extracted from mixed phase. However, as the thermodynamic quantity curve of pure career gas phase is the linear function of temperature, so we may say that such result obtained by subtraction is the pure thermodynamic quantity of Ar.

Figure 7. Assembly value and career gas value *U* obtained by MD method, as a function of temperature.

Figure 8. Assembly value and career gas value *PV* obtained by MD method, as a function of temperature.

Figures 9, 10 compare the results of *U*_{e} and *P*_{e}*V* in MD method obtained based on Figures 7, 8 with such results obtained by using MC method. As already stated in section 2, *U* and *PV* are constituted by two elements, i.e. ideal gas terms and pure interaction terms, see eq. (4).

Upon comparison of MD method and MC method, interaction energy plays an important role. Therefore, the thermodynamic quantities were compared as interaction terms in Figures 9 - 14. The foregoing considered, when we refer to Figures 9, 10, both curves show similar tendency except that *P*_{e}*V* fluctuation on the MD method curve is relatively greater than MC method.

Figure 9. Comparison of *U*_{e} values per one Ar molecule using MD and MC method. (*N*=8)

Figure 10. Comparison of *P*_{e}*V* values per one Ar molecule using MD and MC method. (*N*=8)

In Figure 11, the interaction term of enthalpy *H*_{e} using MD method in accordance with eq. (5) and fitting curve obtained from *H*_{e} are shown.

Figure 11. *H*_{e} value per one Ar molecule using MD method and fittint curve. (*N*=8)

The fitting below is necessary to calculate the interaction term of Gibbs free energy *G*_{e}, because the calculation procedure to obtain *G*_{e} uses fitting function. In addition, Figure 11 does not include the MC method curve because the interaction term of Helmholtz free energy *A*_{e} making pair with *G*_{e} has been calculated for MC method, therefore *H*_{e} is not required to obtain *A*_{e}.

In Figure 12, we obtained two kinds of interaction terms of heat capacities *Cp*_{e} and *Cv*_{e}. The heat capacity at constant pressure *Cp*_{e} is calculated by MD method and another heat capacity at constant volume *Cv*_{e} is calculated by MC method. Each heat capacity was calculated by partial differential in the following eq. (6).

Figure 12. Comparison of heat capacities per one Ar molecule using MD method and MC method.

Although their peaks are not in the same position, such difference is not problematic since they correspond to each phase transition temperature variation.

Interaction term of entropy *S*_{e} by usage of numerical integration described in eq. (7) below is put in Figure 13:

Figure 13. Comparison of *S*_{e} values per one Ar molecule using MD method and MC method. (*N*=8)

In eq. (7), the temperature interval D*T* is constant and it shall be set in accordance with the equation D*T*=*T*_{1}. Accordingly, the first term shall be as (*C*_{e1}/*T*_{1})×*T*_{1}=*C*_{e1}, and *S*_{e} curve and *C*_{e} curve shall have the same starting point.

MD method and MC method were compared using relative value [5] defined by *S*_{e1}= *C*_{e1}. The MD method has a relatively higher *S*_{e} value than the MC method because their values are based on difference of heat capacity *Cp*_{e} and *Cv*_{e}. As the MD method takes into consideration the energy contribution from *P*_{e}*V*, the *S*_{e} value has become higher than that obtained by using MC method.

In Figure 14, comparison of free energy interaction terms obtained from MD method and MC method is indicated.

Figure 14. Comparison of Free Energies per one Ar molecule using MD method and MC method. (*N*=8)

As entropy is used also in eq. (8), the relative values of free energy were compared in the same manner described in the preceding paragraph. There is a difference between constant pressure and constant volume conditions in these two energies. It was able to confirm that phase transition from microcluster to vapor is a voluntary transition because both free energy values show similar curves whose slopes incline in a more negative direction after phase transition. Consequently, we shall conclude that the simulations results of MD method and MC method correspond to each other.

[ 2]

[ 3] Kenji Yasuoka, Mitsuhiro Matsumoto,

[ 4] M.P. Allen, D.J. Tildesley,

[ 5] J.M. Haile,

Return