{'fraction_non_alphanumeric': 0.02410857287364073, 'fraction_numerical': 0.011464216471381607, 'mean_word_length': 4.577809120827457, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Applicability of molecular dynamics method to the prediction of the melting point of refractory metals and compounds

Various methods for measuring the melting point using embedded atom model (EAM) and from first principles modeling are considered. The features of each of them are shown. The results of the work of these methods are compared with the experimental values on the example of systems consisting of pure metals such as Fe, Pd, Ti, Pt Zr, Mo, Ta, W. Good agreement of the estimated melting points with experimental values is shown.


Introduction
Many modern industries, such as the aerospace industry, energetics, transport, need inexpensive refractory materials capable to work in extreme conditions. As known, the most refractory pure metal is tungsten (3695 K), and the most refractory of non-metals is carbon (slightly less than 4000 K). Currently the most refractory compound is Ta 4 HfC 5 , with a melting point of 4215 K [1].
However, empirical methods for researching such materials are complicated by the too large variability of components and the complexity of the experimental measurement of high temperatures. With the increasing computing power in the last few decades, molecular dynamics simulation (MDS) has established itself as a great tool for reasonably accurate prediction of various material properties [1,2]. Thus, MDS is becoming more and more promising to estimate the melting temperature of a material before the start of an experimental study.
In the literature, there are various meanings of the concept of refractoriness, we will understand by this those materials whose melting point is higher than the melting point of iron (1812 K), which is the basis of many alloys.
In this paper, we compared various methods for estimating the melting temperature for pure metals such as Fe, Pd, Ti, Pt, Zr, Mo, Ta, W. The melting temperatures were measured using two different approaches: using the EAM and from first principles (ab-initio). Some simple compounds, such as TaW, TiC and TaC, were also considered.

Methods
One of the distinguishing features of almost any simulation is the idealization of the modeled situation. For example, in molecular dynamics (MD) if we take a box with periodic boundary conditions (in molecular dynamics, this means that if an atom at some integration step flies out one of the edges of the calculation area, then it will be returned to the system from the opposite edge) consisting of tungsten in bcc lattice structure and gradually begin to raise the temperature in it in the iso-enthalpic-isobaric (NPH) ensemble, then it will begin to melt at a temperature significantly higher than it should. This phenomenon has a physical meaning. The fact is that melting requires some center of nucleation. In real life irregularities in the structure or boundary of the sample act as such centers. In this regard, various procedures are used in MD to measure the melting point.

Methods for EAM modeling
For EAM modeling, the LAMMPS software package was used [3]. Interaction potentials were taken from [4]. In all EAM calculations, an integration step of 1 femtosecond was used.
The first approach is to gradually increase the temperature in spheres. As a model system, spheres of different radii (2 nm, 3 nm, 4 nm, 5 nm, 6 nm,) consisting of crystal lattices corresponding to each metal at 300 K were considered. The temperature of the spheres increased at a rate of 10 12 K/s, before becoming a completely liquid drop, as shown in figure 1 for a TaW particle. Figure 1a shows the initial structure of a TaW particle at 300 K. Figure 1b shows the process of melting a particle (at a temperature of 3574 K). It is clearly seen that melting goes from the boundary of the particle to the center. Figure 1c shows a fully melted drop of TaW. The melting point for each particle was calculated as the temperature at which the enthalpy jump occurred in the system. As is known, particles with a small radius will have a lower melting point. To calculate the melting point of the bulk metal, linear approximation of the values obtained for spheres of different radii to a sphere of infinite radius was carried out (1 / (Radius of the sphere) values were taken along the X axis and enthalpy jump temperatures measured along the Y axis). For example, for TaW particles with radii of 2 nm, 4 nm, 6 nm, we have 3000 K, 3450 K, 3574 K, respectively. Approximation to a sphere of infinite radius gives the melting temperature for TaW equal to 3872 K, which is 300 K higher than the experimental value. The melting temperature T m eam for pure metals obtained by this method is shown in table 1.
The second approach is known as the solid-liquid coexistence method. The idea of the method consists in the equilibrium coexistence of the solid and liquid phases of a substance in one simulated system at a certain temperature. For this, a simulated box of 50 nm in x direction, 10 nm in y direction and 5 nm in z direction was created. Periodic boundary conditions along axes y and z were used. Initially, the lattice corresponding to the compound was set in the simulated system. Then the system was divided into the left and right equally large half. The left half of the sample was heated to a temperature above the melting point in the isothermal-isobaric ensemble (NPT) until completely liquid phase is reached, while the right half remained untouched in the iso-enthalpic-isobaric (NPH) ensemble at room temperature. A typical model system is shown in figure 2. Then we need to guess the approximate melting point for a metal or compound, let's call this value T i .  Figure 2. The structure of the coexistence of solid and liquid phases on the example of TaW system. The blue spheres correspond to tantalum atoms. The yellow spheres correspond to tungsten.
In each of the halves, we set the temperature T i in the NPT ensemble. Then the simulation is again carried out in the whole system in the NPH ensemble. Depending on how closely we guessed the T i , there may be different behavior in the evolution of such a system. If the phase boundary remains stationary, the melting point is equals to T i . If the phase boundary shifts towards the crystalline component, but at the same time there remains at least some crystalline phase in the system, the melting point will be equal to the new steady-state temperature in the system. Similarly, if part of the liquid crystallizes at the phase boundary, then the melting point will be equal to the new steady-state temperature in the system. It is also possible that the whole system become completely liquid or a completely crystallized. In these cases, we can only say that the melting point is lower or higher, respectively, than the chosen initial temperature. And the simulation needs to be done again with the modified value of T i . With this approach, the larger the size of the simulated system, the less accurately you need to be to guess the initial temperature T i . The values obtained in this way agree fairly accurately with the values T m eam .

Methods for ab-initio modeling
For ab-initio modeling, the VASP software package was used. Model systems contained from 100 to 200 atoms. PAW (projector augmented wave) pseudopotentials were used for the simulation. The integration step was equal to 1 femtosecond.
The methods used for EAM modeling are only partially suitable for first-principle modeling. This is due to the fact that computational complexity in ab-initio modeling grows as the third power of the number of valence electrons in the system. Therefore, long-time system modeling with more then 500 atoms using first-principle methods is still considered challenging task even for powerful computing systems.
The first approach to measure the melting point is to gradually heat the system in canonical ensemble (NVT) until reaching a completely liquid state. For this for each metal, a crystal lattice corresponding to it was set at 300 K and simulation was performed for 0.1 picoseconds and then the temperature increased abruptly by 50 K and calculations were again performed for 0.1 picoseconds and so continued until a fully molten state was reached. The fact that the volume remains constant makes the melting process easier, however, the melting point obtained in this way will correspond to a rather strongly increased external pressure. Such an approach can only serve as a rough estimate of the melting point.
The second approach for a more accurate calculation of the melting point is a variation of the solidliquid coexistence method. As well as during EAM modeling, a system consisting of a certain amount of a substance in a liquid and crystalline states was created. Due to computational limitations, the number of atoms in the system was less than 200. For such a small system, reaching an equilibrium state in which both the crystalline and liquid phases would remain is almost impossible. Thus for each compound, a series of experiments was carried out in some of which, the system went into a completely liquid state, and the others became completely crystalline. Temperature equilibrium of a liquid and a solid with such an approach will be that temperature at which the transition to both completely liquid and fully crystalline state is statistically equally likely. The melting points obtained by this method for pure metals T m ab-in are listed in the table 1. It should be taken into account that the error of this approach is about 50 K, but can be reduced by increasing the statistical base, i.e., carrying out a larger number of computational experiments.
Whether the structure is liquid or crystalline was determined by the method of CNA (common neighbor analysis) [5], which allows to determine the belonging of atoms to a cubic face-centered lattice, a cubic body-centered lattice and a hexagonal close-packed lattice (and all the pure metals we consider have one of the above structures). The system was considered liquid when the number of atoms defined as crystalline by CNA method became less than 5%.

Results and discussion
The experimental melting points T m * are taken from [6]. As can be seen, the use of EAM potentials gives rather far melting points from real values. The difference between T m eam and T m * is from 74 K to 443 K. As can be seen, the difference between the temperature T m ab-in and the experimental values is about 50 K, which fits into the error of the calculation method. The EAM modelling requires potentials of cross-pair interaction that can not always be found in the literature. For example, the calculated melting point for titanium carbide using the potential from work [7] was equal to 3400 K, which is consistent with the experimental data. And the melting point of tantalum carbide calculated by the potential developed by a method similar to the method in work [7] was equal to 3800 K, which is about 350 degrees lower than the experimental value. Ab-initio modeling using the method described above gives a melting point 3450 K for TiC and 4150 for TaC, which agrees with experimental values, with an error of 50 degrees.

Conclusions
EAM modeling is widely used to study the processes occurring in some systems, including it can be useful for studying the processes of synthesis and consolidation of materials. This is due to the fact that it is possible to model systems consisting of millions of atoms in a reasonable time. But for accurate measurement of the melting point, like many other physical quantities, EAM modeling is not the best choice. Ab-initio simulation gives melting points quite close to experimental data. The main disadvantage is the smallness of the simulated system, due to the increasing complexity of the calculation. Another advantage of ab-initio modeling is the presence of pseudopotentials for all elements of the periodic system.