Hugoniot States and Mie–Grüneisen Equation of State of Iron Estimated Using Molecular Dynamics

: The objective of this study was to develop a micromechanical approach for determining the Mie–Grüneisen EOS parameters of iron under the Hugoniot states. The multiscale shock technique (MSST) coupled with molecular dynamics (MD) simulations was employed to describe the shocked Hugoniot relation of single-crystal (SC) and nanocrystalline (NC) iron under high pressures. The Mie–Grüneisen equation of state (EOS) parameters, the cold pressure ( P c ), the cold energy ( E c ), the Grüneisen coefﬁcient ( γ ), and the melting temperature ( T m ) are discussed. The error between SC and NC iron results was found to be less than 1.5%. Interestingly, the differences in Hugoniot state ( P H ) and the internal energy between SC and NC iron were insigniﬁcant, which shows that the effect of grain size (GS) under high pressures was not signiﬁcant. The P c and E c of SC and NC iron calculated based on the Morse potential were almost the same with those calculated based on the Born–Mayer potential; however, those calculated based on the Born–Mayer potential were a little larger at high pressures. In addition, several empirical and theoretical models were compared for the calculation of γ and T m . The Mie–Grüneisen EOSs were shown on the 3D contour space; the pressure obtained with the Hugoniot curves as the reference was larger than that obtained with the cold curves as the reference.


Introduction
Iron is a major element in the Earth's core. The properties of iron at elevated temperatures and pressures are crucial to model the temperature profile and the dynamics of the earth core, where pressures could reach up to 130-370 GPa [1,2]. Under high-pressure conditions, the behavior of metallic materials can be described by the equation of state (EOS) [3]. In the study on high-pressure physics, a deep understanding of EOS and other characteristics, such as the shock Hugoniot, Grüneisen coefficient (γ), and melting temperature (T m ), is vital [4,5].
The shock Hugoniot relations, which can be applied to calculate the EOS, have received great attention [6,7]. Barker [8] reexamined the α-iron Hugoniot at around 13 GPa. The results showed that the Hugoniot α-phase iron was quite normal under 13 GPa, and there was no indication of an early phase transition to -iron. Thomas et al. [9] investigated the shock Hugoniot of α-iron under 10 GPa through low-velocity symmetric impact experiments. The results provided reliable data and fit for α-iron under low-stress conditions below 10 GPa. Above 13 GPa, a solid-solid phase (α-) transition will occur in iron [10]. Brown et al. [11] collected the shock Hugoniot data of iron (<0.3% carbon) in Los Alamos 2 of 17 from about 39 to 442 GPa, which were obtained using explosive driven shock and a twostage light-gas gun system. Zaretsky [12] studied the shock response of polycrystalline iron using the planar impact experiment method and determined the borders of stability of iron αand γ-phases. Liu et al. [13] investigated the crystal orientation effect on the phase transition stress of single iron under shock compression; the results show that the orientation effect on the phase transition stresses of single crystal iron is not obvious. Sjostrom and Crockett [14] conducted quantum molecular dynamics simulation to investigate the shock Hugoniot and equation of state of warm dense iron; the shock Hugoniot of up to 20 TPa agrees with the experimental results. Kadau et al. [15,16] investigated the shock-induced phase transformations in single-crystal (SC) and polycrystalline iron using nonequilibrium molecular dynamics (NEMD) simulations. The dynamics of shock-induced transformations were discussed, and the Hugoniot data of SC and polycrystalline iron under pressures up to 300 GPa were obtained.
The Mie-Grüneisen EOS, which relates the pressure (P), specific volume (V), and internal energy (E), is widely used to characterize the behavior of metal materials in highpressure physics [17][18][19][20][21][22][23]. In the study by Heuzé [24], the general P(V, E) Mie-Grüneisen form was extended to a complete S(V, E) form, which could provide all of the thermodynamic properties. Zhang et al. [25] proposed a cold energy mixture theory in which the Mie-Grüneisen EOS is calculated from the cold curves and Hugoniot curves, and the authors also examined γ. In the Mie-Grüneisen EOS, γ plays an important role, as it can clarify the melting of solids at extreme conditions [26]. There are three models for characterizing the Grüneisen coefficient γ, the Slater model [27], the Dugdale-MacDonald model [28], and the free-volume model [29]. In early 1963, Vashchenko and Zubarev [29] proposed the idea of writing a unified equation for the three models. Based on this, Al'tshuler [30] gives a unified expression of the Grüneisen coefficient γ. Cui and Yu [31] proposed a model to calculate γ at high temperatures. The model agrees well with the experimental data of iron in a wide range of pressures. Sha and Cohen [32] used first-principle calculations to study the γ of bcc iron under pressures up to 40 GPa and temperature of 2000 K. The calculated γ showed good agreement with the experiments, and the variation of γ with temperature was remarkable, but the variation with pressure was moderate. Li et al. [33] calculated the γ of iron in the pressure range of 90-160 GPa based on the Hugoniot data obtained by plate impact experiments. The results show that γ 0 ·ρ 0 = γ·ρ = const in the pressure range.
Shock-induced melting is an important phenomenon for studying the solid-liquid phase transition and the nature of material melting at extreme conditions, and the melting temperature T m is related to γ 0 . Based on the Lindemann melting criterion and the volume dependence of the Gruneisen coefficient γ, Yang et al. [5] investigated the melting characteristics of single and polycrystalline Al. Errandonea [34] used a laser-heated diamond-anvil cell to study the melting characteristics of Cu, Ni, Pd, and Pt in high-pressure conditions; it was shown that Simon equations could be used to describe the melting curves. Anzellini et al. [35] conducted static laser-heated diamond anvil cell experiments to investigate the melting point of iron under pressures up to 200 GPa. The results indicated that the temperature of iron at the inner core boundary was 6230 ± 500 K. Bouchet et al. [36] calculated the melting curve, and EOS of hcp and bcc iron under pressures up to 1500 GPa; the melting temperature was 11,000 K at the highest pressure.
MD simulation has become an important research method in high-pressure physics. For example, the NEMD method [37], equilibrium molecular dynamics (EMD) method, [38] and first-principle method [39] have been successfully applied to study the shock Hugoniot, γ, and T m . The time scale and model size in NEMD and EMD simulations are restricted by computing resources; as the time scale and model size increase, the computational time will significantly increase. The first-principle method can fundamentally calculate the molecular structure and material properties, but the model size is too small in the shockwave propagation simulation. To solve the computational time and model size shortcomings in the above methods, Reed et al. [40] proposed a multiscale simulation method, called the multiscale shock technique (MSST), to study the shockwave propagation in materials. Compared with the above methods, the MSST method can efficiently save computational time and guarantee simulation accuracy with an acceptable model size of the MD system [41,42]. Thus, in the present study, the MSST method is applied to obtain the shock Hugoniot data in SC and nanocrystalline (NC) iron.
In the literature [3,9,12,13,15,16], there are relatively few studies on the comparison between the shock Hugoniots, and EOS of SC and polycrystalline iron. The EOS, shock Hugoniot, and T m are significant physical properties that must be examined to study the behavior of metal materials under high-pressure conditions. There are two issues involved in the relation between shock Hugoniots and EOS: one regards obtaining the EOS based on the shock Hugoniots, and the other concerns the determination of the shock Hugoniots based on the EOS. In addition, the Grüneisen parameter and melting temperature can be calculated from the shock Hugoniot data. The investigation path of this work firstly utilizes MD simulation to obtain the shock Hugoniots of SC and NC iron; thus, the Hugoniot curve can be obtained. Therefore, based on the Hugoniot data, γ and T m can be determined. Finally, the EOS of SC and NC iron can be presented.

Hugoniot Pressure (P H ) and Internal Energy (E H )
The Hugoniot relation can be expressed by three conservation equations: the mass equation, the momentum equation, and the energy equation: where ρ, u s , u p , P H , E H , and V (=1/ρ) are the density, shockwave velocity, particle velocity, pressure, and internal energy per unit mass and specific volume of the shocked material, respectively. The shock wave velocity in the MD simulation was set to 6-11 km/s, and the particle velocity under different loading conditions was obtained by varying the loading velocity. The subscript "0" denotes the qualities in the initial state without shock, and "H" represents the quantities in the Hugoniot state. For example, P H represents the pressure in the Hugoniot state. It has been found that under a wide range of pressure, the relationship between the shockwave velocity u s and the particle velocity u p is approximately linear [11]: where C 0 and λ are the volume sound speed at zero pressure and the fitting parameter, respectively. From MD simulation results, C 0 and λ of SC and NC iron can be determined. Once C 0 and λ is determined, it can be used to calculate the Hugoniot curve P H (V) and the internal energy E H (V). According to Equations (1)-(4), when P 0 = 0, the Hugoniot curve P H (V) and the internal energy E H (V) can be derived: where c v and T 0 are the specific heat at constant volume and the initial temperature, respectively. If the E H is much larger than E 0 , E 0 can usually be ignored.

Cold Pressure (P c ) and Cold Energy (E c )
From the Hugoniot relations at room temperature (usually 300 K), V 0K , ρ 0K , C 0 , and λ can be obtained as follows [43]: where γ(V 0 ) is the Grüneisen coefficient at room temperature (usually 300 K); α v (T) is the volume expansion coefficient, which is approximately constant from 0 K to 300 K; and C 0 , and λ are parameters of u s and u p relations at 0 K. The Born-Mayer [44] and Morse [45] potentials are usually applied to describe the atomic forces at 0 K; thus, the cold pressure P c and cold energy E c can be expressed as follows: where Q, q, A, and B are material constants; δ = V 0K /V; and ρ 0K is the material density at 0. At P = 0 and T = 0, we can assume that the first-order and second-order derivatives of P H , isentropic P s , and P c are almost equal [5,46]. Thus, we have: With the determined C 0 and λ from MD simulations, and according to Equations (5) and (10)-(15), Q, q, A, and B can be calculated as follows:

Grüneisen Coefficient γ
The Grüneisen coefficient γ shows a dependence of volume and can be expressed as follows [30]: where t is used to represent three models for characterizing γ. When t = 0, Equation (23) represents the Slater model [27] γ s (V); t = 1 represents the Dugdale-MacDonald model [28] γ DM (V); and t = 2 represents the free-volume model [29] γ f (V). At P = 0 and V = V 0 , Equation (22) approximately satisfies the isentropic P s (V), and, consequently, Equation (22) can be expressed as follows: When t is equal to 0, 1, and 2: When P = 0 and V = V 0 , the first-order and second-order derivatives of P s (V) and P H (V) are identical. In addition, P H (V 0 ) and P H (V 0 ) can be derived from P H obtained via MD simulations or experiments; therefore, we have: Substituting Equations (27)- (29) into Equation (23) yields As we have determined λ from MD simulations, the value of Equations (30)-(32) can be obtained. From previous work, the volume dependence of γ can be commonly expressed by the following simpler empirical formulas [5,[47][48][49]:

Melting Temperature (T m )
In high-pressure physics, the melting temperature T m (V) is a significant material characteristic. The Lindemann melting criterion is given by [50]: According to Equation (39), Equations (40)-(44) are deduced as: where γ 0 is determined from Equations (30)-(32), T m0 is the melting temperature at 0 GPa, and is 1811 K for iron [51].

Mie-Grüneisen Equation of State
The Mie-Grüneisen EOS is usually applied to describe the behavior of metal materials under extreme shock compression, and the form of Mie-Grüneisen EOS can be expressed as follows [24]: From the calculation results of P H , E H , P c , E c , and γ(V), the Mie-Grüneisen EOS could be obtained.

MD Simulation
The embedded-atom method (EAM) potential presented by Mendelev et al. [52] was used to conduct the MD simulation in this study. Three simulation models were established to investigate the dynamic characteristics of SC and NC iron under intensive shock conditions. The established models were based on bcc iron, and the shockwave direction is [001] in SC iron. Figure 1 shows the three simulation models. The NC models were built on the Voronoi algorithm in Atomsk [53]. The periodical boundary conditions were set in all directions for both SC and NC iron. To balance the calculation time and model size, the simulation system of SC iron was about 8.4 × 8.4 × 8.4 nm 3 with 54,000 atoms. The model sizes of NC iron were about 10.9 × 10.9 × 10.9 nm 3 with 111,195 atoms and 54.7 × 54.7 × 54.7 nm 3 with 8,268,098 atoms. There were two grain sizes (GSs): 5 nm and 25 nm, which were chosen to study the GS effect in the shocking behavior.
conditions. The established models were based on bcc iron, and the shockwave direction is [001] in SC iron. Figure 1 shows the three simulation models. The NC models were built on the Voronoi algorithm in Atomsk [53]. The periodical boundary conditions were set in all directions for both SC and NC iron. To balance the calculation time and model size, the simulation system of SC iron was about 8.4 × 8.4 × 8.4 nm 3 with 54,000 atoms. The model sizes of NC iron were about 10.9 × 10.9 × 10.9 nm 3 with 111,195 atoms and 54.7 × 54.7 × 54.7 nm 3 with 8,268,098 atoms. There were two grain sizes (GSs): 5 nm and 25 nm, which were chosen to study the GS effect in the shocking behavior.
First, the simulation system was equilibrated at 300 K for 10 ps in the NPT ensemble to ensure that the system reached a steady state in order to eliminate the residual stress. Then, the shockwave was examined using the MSST method under velocities of 6-11 km/s. Moreover, the computation time should be long enough to guarantee that the shockwave in the system reaches a steady state after the shock. Thus, the statistical results, such as the particle velocity, pressure, and temperature, could be obtained from the MD results.

Shock Hugoniot Pressure and Internal Energy
According to the experimental Hugoniot data from Brown et al. [11], the shock velocity us ranged from 5.3 to 11.5 km/s at about 40-440 GPa. To ensure the stability of the MD simulation and evaluate the model validity, the velocity of us was set in the range of 6 to 11 km/s in this study. Figure 2 shows the relationship between us and up for SC iron obtained via the MSST method; there is a good agreement between the MD simulation results and experimental data. Moreover, the linear relationship between us and up is in excellent accordance with the data in Figure 2. The linear relationship is also consistent with the study by Prieto and Renero [54], which indicates that the relationship between us and up of iron below 5 megabars can be described by a linear relationship. For SC iron, the relationship between us and up is given as us = 4.071 + 1.538up. The relationship between us and up from the experimental data is expressed as us = 3.935 + 1.578up. The relative errors First, the simulation system was equilibrated at 300 K for 10 ps in the NPT ensemble to ensure that the system reached a steady state in order to eliminate the residual stress. Then, the shockwave was examined using the MSST method under velocities of 6-11 km/s. Moreover, the computation time should be long enough to guarantee that the shockwave in the system reaches a steady state after the shock. Thus, the statistical results, such as the particle velocity, pressure, and temperature, could be obtained from the MD results.

Shock Hugoniot Pressure and Internal Energy
According to the experimental Hugoniot data from Brown et al. [11], the shock velocity u s ranged from 5.3 to 11.5 km/s at about 40-440 GPa. To ensure the stability of the MD simulation and evaluate the model validity, the velocity of u s was set in the range of 6 to 11 km/s in this study. Figure 2 shows the relationship between u s and u p for SC iron obtained via the MSST method; there is a good agreement between the MD simulation results and experimental data. Moreover, the linear relationship between u s and u p is in excellent accordance with the data in Figure 2. The linear relationship is also consistent with the study by Prieto and Renero [54], which indicates that the relationship between u s and u p of iron below 5 megabars can be described by a linear relationship. For SC iron, the relationship between u s and u p is given as u s = 4.071 + 1.538u p . The relationship between u s and u p from the experimental data is expressed as u s = 3.935 + 1.578u p . The relative errors of C 0 and λ obtained from the MSST method and experiment were 3.3% and 2.6%, respectively, which proves the accuracy of the MD simulation.
To study the GS effect on the u s -u p relationship for NC iron, MD simulation was conducted at the same shock velocity as that of SC iron. The results are shown in Figure 3. The results of NC iron with different GS values are close. These discrete data are characterized by u s = 4.011 + 1.552u p and u s = 4.01 + 1.553u p for GS 5 nm and GS 25 nm, respectively. The relative error of C 0 and λ between GS 5 nm and GS 25 nm was less than 0.1%, and compared with the SC iron data, C 0 and λ were approximately equal to an error of <1.5%. The results indicate that the GS effect had little influence on the linear relationship between u s and u p for SC and NC iron, and it could be neglected in the case of NC iron. of C0 and λ obtained from the MSST method and experiment were 3.3% and 2.6%, respectively, which proves the accuracy of the MD simulation. To study the GS effect on the us-up relationship for NC iron, MD simulation was conducted at the same shock velocity as that of SC iron. The results are shown in Figure 3. The results of NC iron with different GS values are close. These discrete data are characterized by us = 4.011 + 1.552up and us = 4.01 + 1.553up for GS 5 nm and GS 25 nm, respectively. The relative error of C0 and λ between GS 5 nm and GS 25 nm was less than 0.1%, and compared with the SC iron data, C0 and λ were approximately equal to an error of <1.5%. The results indicate that the GS effect had little influence on the linear relationship between us and up for SC and NC iron, and it could be neglected in the case of NC iron. The shock Hugoniot data obtained from MD simulations of SC and NC iron are listed in Table 1. As can be seen in Table 1, the particle velocity is high in NC iron below a shock velocity of 8.0 km/s, with a maximum difference of 2.5% between SC and NC iron. However, as the shock speed continues to grow, the difference in particle velocity becomes  To study the GS effect on the us-up relationship for NC iron, MD simulation was conducted at the same shock velocity as that of SC iron. The results are shown in Figure 3. The results of NC iron with different GS values are close. These discrete data are characterized by us = 4.011 + 1.552up and us = 4.01 + 1.553up for GS 5 nm and GS 25 nm, respectively. The relative error of C0 and λ between GS 5 nm and GS 25 nm was less than 0.1%, and compared with the SC iron data, C0 and λ were approximately equal to an error of <1.5%. The results indicate that the GS effect had little influence on the linear relationship between us and up for SC and NC iron, and it could be neglected in the case of NC iron. The shock Hugoniot data obtained from MD simulations of SC and NC iron are listed in Table 1. As can be seen in Table 1, the particle velocity is high in NC iron below a shock velocity of 8.0 km/s, with a maximum difference of 2.5% between SC and NC iron. However, as the shock speed continues to grow, the difference in particle velocity becomes  The shock Hugoniot data obtained from MD simulations of SC and NC iron are listed in Table 1. As can be seen in Table 1, the particle velocity is high in NC iron below a shock velocity of 8.0 km/s, with a maximum difference of 2.5% between SC and NC iron. However, as the shock speed continues to grow, the difference in particle velocity becomes very small. At 11 km/s, the maximum shock pressure in SC and NC iron is 387.98-389.62 GPa. From [11], it can be found that the shock pressure is 361.8 GPa at a shock velocity of 10.7 km/s and 400.7 GPa at a shock velocity of 11.11 km/s. The MD results are within the range of the experimental data, which shows the accuracy of the simulation.
A detailed comparison of MD simulation results and experimental data is presented in Table 2. From Table 2, it can be observed that the volume sound speed C 0 in SC iron is higher than that in NC iron; however, the fitting parameter λ in SC iron is lower than that in NC iron, and the same characteristics are also observed in SC and NC Al [5]. The maximum relative error of C 0 between the experiment and SC iron is 3.34%, and the maximum relative error of NC iron is 1.89%. The difference of C 0 between SC and NC iron is~1.5%, which is less than the difference with experiment data. The same characteristics in λ can be observed; the maximum relative error between the experiment and SC iron is 2.60%, and the maximum relative error of NC iron is 1.67%. The difference of λ between SC and NC iron is~0.96%. The results indicate that a shock response difference between SC and NC iron does exist, but the difference is not very significant.  Table 2. C 0 and λ obtained from the linear relationship between u s and u p for SC and NC iron at 300 K. The initial density is 7.85 g/cm 3 in both the experiment and MD simulation.

Material Shock Velocity
C 0 (km/s) λ Based on the linear relationship obtained for SC and NC iron, the Hugoniot curve could be obtained using Equations (5) and (6). Figure 4 shows the Hugoniot curves of SC and NC iron. The black circle represents the results by Brown et al. [11], and the colored curves represent this work for SC and NC iron. Figure 4a shows that for V/V 0 > 0.6, the curves obtained through the MSST method agree well with the results of Brown et al. [11]; when V/V 0 < 0.6, the curves are approximately equal to the experiment results, and the deviation is less than 3%. Considering the insignificant error of C 0 and λ for SC and NC iron, there is no noticeable difference between the three curves in Figure 4a. As shown in Figure 4b, the internal energies obtained in the MD simulation and experiment were almost the same; when V 0 /V = 1.7, the internal energy of SC iron was smaller than that of NC iron, and the deviation at V 0 /V = 1.7 was 0.5%. The results not only prove the reliability of the linear relationship between u s and u p but also show the EAM potential for iron in the MD simulation. In addition, the results signify that the MD results are appropriate for describing the cold pressure, the cold energy, the melting temperature, and the Mie-Grüneisen EOS. Figure 4b, the internal energies obtained in the MD simulation and experiment were almost the same; when V0/V = 1.7, the internal energy of SC iron was smaller than that of NC iron, and the deviation at V0/V = 1.7 was 0.5%. The results not only prove the reliability of the linear relationship between us and up but also show the EAM potential for iron in the MD simulation. In addition, the results signify that the MD results are appropriate for describing the cold pressure, the cold energy, the melting temperature, and the Mie-Grüneisen EOS.

Cold Pressure and Cold Energy
Based on the C0 and λ obtained from the Hugoniot data, the material constants for calculating cold pressure and cold energy in Equations (18)-(21) could be obtained. The values are presented in Table 3.

Cold Pressure and Cold Energy
Based on the C 0 and λ obtained from the Hugoniot data, the material constants for calculating cold pressure and cold energy in Equations (18)-(21) could be obtained. The values are presented in Table 3.  Figures 5 and 6 compare the cold pressure P c and cold energy E c between SC and NC iron, with the experimental data as a reference. The solid and dashed curves represent P c and E c obtained using the Born-Mayer and Morse potentials, respectively. The results from the Born-Mayer and Morse potentials presented the same increasing trend: P c and E c increased with volume compression (V 0K /V), which presented a sharp increment at higher compression. of the potential function used, the cold pressure curve of SC iron was above that of the NC iron. This indicates that the cold pressure of SC iron was higher than that of NC iron, although the difference was quite small. Under V0K/V of less than 1.4, the cold energy curves agree well with each other, as the cold pressures obtained by the Born-Mayer and Morse potentials in the low-compression region are approximately equal. Under V0K/V of above 1.4, the Ec curves gradually deviate from each other. The cold energy obtained by the Born-Mayer potential is larger than that from the Morse potential, and the same trend was found for the cold pressure. The relative error of Ec is about 7% when V0K/V is 2.0. Furthermore, the cold energy in the SC iron is larger than that of the NC iron in both Born-Mayer and Morse potentials.

Grüneisen Coefficient
The material constant λ in SC and NC iron was used to calculate the γ0 of the Slater, Dugdale-Macdonald, and free-volume models through Equations (30)- (32). Subsequently, γ0(V) was substituted into Equations (34)- (38) to obtain γ(V), and the experimental data were used for comparison. Figure 7 shows the results of γ(V) from different empirical and theoretical models. At V 0K /V less than 1.3, the curves of P c agree well with each other. The result shows that the cold energies described by the Born-Mayer and Morse potentials are basically the same in the low-compression region. When the V 0K /V is above 1.3, the P c curves gradually deviate from each other, and the P c curves obtained by the Born-Mayer potential gradually become larger. The relative error of P c was about 10% when V 0K /V was 2.0. Irrespective of the potential function used, the cold pressure curve of SC iron was above that of the NC iron. This indicates that the cold pressure of SC iron was higher than that of NC iron, although the difference was quite small.
Under V 0K /V of less than 1.4, the cold energy curves agree well with each other, as the cold pressures obtained by the Born-Mayer and Morse potentials in the low-compression region are approximately equal. Under V 0K /V of above 1.4, the E c curves gradually deviate from each other. The cold energy obtained by the Born-Mayer potential is larger than that from the Morse potential, and the same trend was found for the cold pressure. The relative error of E c is about 7% when V 0K /V is 2.0. Furthermore, the cold energy in the SC iron is larger than that of the NC iron in both Born-Mayer and Morse potentials.

Grüneisen Coefficient
The material constant λ in SC and NC iron was used to calculate the γ 0 of the Slater, Dugdale-Macdonald, and free-volume models through Equations (30)- (32). Subsequently, γ 0 (V) was substituted into Equations (34)- (38) to obtain γ(V), and the experimental data were used for comparison. Figure 7 shows the results of γ(V) from different empirical and theoretical models.
Equations (34)- (38) are commonly used empirical models for calculating γ(V) at high pressure. From Figure 7, it can be seen that the Grüneisen coefficient γ(V) first rapidly decreased and then slowly dropped as the V 0 /V increased. At high compression, the difference between SC and NC results became smaller; this indicates that the microstructure has little effect on γ(V) under high pressure. Equation (38) is suitable for low-compression data, that is, V 0 /V up to 1.5, and is free with λ; thus, different models do not affect their results. With Equation (38) as a reference, when V 0 /V = 1.0, the initial value of γ 0 is γ f (V 0 ) < γ DM (V 0 ) < γ S (V 0 ). For all three models, the γ(V) obtained by the experimental data (solid curve) is at the top followed by that for the NC iron (dash curve), and that for the SC iron (dot curve) is at the bottom due to the difference between the λ values. Equations (34)- (38) are commonly used empirical models for calculating γ(V) at high pressure. From Figure 7, it can be seen that the Grüneisen coefficient γ(V) first rapidly decreased and then slowly dropped as the V0/V increased. At high compression, the difference between SC and NC results became smaller; this indicates that the microstructure has little effect on γ(V) under high pressure. Equation (38)   Equation (38) (magenta curve) is free of λ, and it has been proven to agree well with the experimental data under V 0 /V up to 1.5 [49]. The γ(V) obtained using the free-volume model has the smallest difference under V 0 /V up to 1.5. When V 0 /V is 1.5, the γ(V) obtained by Equation (37) (cyan curve) intersects with the magenta curve. This indicates that Equation (27) is applicable to characterize the γ(V) of iron at higher compression compared with other expressions. This has also been proved by Jacobs and Schmid-Fetzer [48]. In Figure 7c, the green curve is closer to the cyan curve than other curves. The result shows that Equation (36) (green curve) is more suitable than other expressions for describing γ(V). The relative error of the green curve is 7.8% at V 0 /V = 2.0.
< γDM(V0) < γS(V0). For all three models, the γ(V) obtained by the experimental data (solid curve) is at the top followed by that for the NC iron (dash curve), and that for the SC iron (dot curve) is at the bottom due to the difference between the λ values. Compared with the results from Equations (35)- (38), the results from Equation (34) (black curves) show a more significant downward trend with the increase in V0/V. This shows that the γ(V) calculated using Equation (34) has the maximum change rate with the increase in V0/V. The results from Equation (35)-(38) present a similar trend in predicting γ(V).
Equation (38) (magenta curve) is free of λ, and it has been proven to agree well with the experimental data under V0/V up to 1.5 [49]. The γ(V) obtained using the free-volume model has the smallest difference under V0/V up to 1.5. When V0/V is 1.5, the γ(V) obtained by Equation (37) (cyan curve) intersects with the magenta curve. This indicates that Equation (27) is applicable to characterize the γ(V) of iron at higher compression compared with other expressions. This has also been proved by Jacobs and Schmid-Fetzer [48]. In Figure 7c, the green curve is closer to the cyan curve than other curves. The result shows that Equation (36) (green curve) is more suitable than other expressions for describing γ(V). The relative error of the green curve is 7.8% at V0/V = 2.0.

Melting Temperature
Based on the Grüneisen coefficient and the Lindemann law, the melting temperature Tm can be obtained. In this study, γ0 is calculated via the free-volume model with Equation (26), and Tm is predicted by Equations (40)- (44). The relationship between Tm and pressure is illustrated in Figure 8, together with experimental data [35,[55][56][57][58][59][60][61][62][63] and theoretical data [64][65][66][67][68]. As shown in Figure 8, the temperature rises quickly to the melting point because of the adiabatic process under shock compression, and Tm shows an increasing trend as the pressure increases. The results from Equations (41)- (42) and (44) present a similar and close trend; the red curve at the top corresponds to Equation (41), while the green curve at the bottom corresponds to Equation (43), and the maximum error is 29.5% at 450 GPa. Moreover, the melting temperature in SC iron is a little lower than that in NC iron; the maximum error is 2.5% at 450 GPa.
At lower pressure (<100 GPa), all curves except the green one are very close; this shows that the prediction abilities of Equations (40)- (42) and (44)   As shown in Figure 8, the temperature rises quickly to the melting point because of the adiabatic process under shock compression, and T m shows an increasing trend as the pressure increases. The results from Equations (41), (42) and (44) present a similar and close trend; the red curve at the top corresponds to Equation (41), while the green curve at the bottom corresponds to Equation (43), and the maximum error is 29.5% at 450 GPa. Moreover, the melting temperature in SC iron is a little lower than that in NC iron; the maximum error is 2.5% at 450 GPa.
At lower pressure (<100 GPa), all curves except the green one are very close; this shows that the prediction abilities of Equations (40)- (42) and (44) in the low-compression area are almost the same. It can be found that T m calculated by Equations (40) and (44) agrees well with DAC experimental data in [57], and T m calculated by Equation (43) agrees well with DAC experimental data in [59,60,62]. T m calculated by Equations (41) and (42) is little higher than that obtained from DFT calculations [64][65][66]. However, the difference between the results gradually increased at higher pressure (>100 GPa). The increasing trend indicates that in the low-compression area, the effect of γ 0 on T m is insignificant; the effect of γ 0 on T m becomes more significant with increasing pressure. This feature shows the opposite trend compared with the γ(V) results. In Figure 8, it can be found that the DAC experimental data of Anzellini et al. [35], the theoretical data obtained by MD calculations conducted by Zhang et al. [58], and DFT calculations conducted by Alfè et al. [64][65][66] are in good agreement. Thus, the orange melting curve obtained by the DAC experiment up to~350 GPa was chosen as the main reference in this study. The results obtained from Equations (40)- (44) are similar to the orange line, although there are some differences in the results. At lower pressure (<100 GPa), the black line obtained by Equation (40) is the closest to the orange line, and the green line (Equation (43)) and purple line (Equation (44)) slightly deviate from the orange line. This indicates that Equation (40) can be used to approximately predict T m under 100 GPa. In the entire pressure range, the orange line is within the calculation result from Equation (44) at pressure less than 300 GPa; the result obtained from Equation (44) is lower than the orange line above 300 GPa. This result shows that the calculation results under 300 GPa can cover the experimental data, and the results obtained at higher pressures are lower than the experimental data.

Mie-Grüneisen Equation of State
The Mie-Grüneisen EOS can be obtained from Equations (45) and (46) with the Hugoniot curve (Equation (45)) and the cold curve (Equation (46)) as the reference curves. The Hugoniot curve is obtained via MD simulation. Based on previous analysis, the Grüneisen coefficient γ in this study is calculated by Equations (32) and (36). The cold pressure and cold energy are calculated from the Born-Mayer potential. Figure 9 shows the Mie-Grüneisen EOS of iron in 3D contour form of SC iron based on the cold pressure and cold energy, with the Hugoniot relations as the reference.
lations conducted by Zhang et al. [58], and DFT calculations conducted by Alfè et al. [64][65][66] are in good agreement. Thus, the orange melting curve obtained by the DAC experiment up to ~350 GPa was chosen as the main reference in this study. The results obtained from Equations (40)- (44) are similar to the orange line, although there are some differences in the results. At lower pressure (<100 GPa), the black line obtained by Equation (40) is the closest to the orange line, and the green line (Equation (43)) and purple line (Equation (44)) slightly deviate from the orange line. This indicates that Equation (40) can be used to approximately predict Tm under 100 GPa. In the entire pressure range, the orange line is within the calculation result from Equation (44) at pressure less than 300 GPa; the result obtained from Equation (44) is lower than the orange line above 300 GPa. This result shows that the calculation results under 300 GPa can cover the experimental data, and the results obtained at higher pressures are lower than the experimental data.

Mie-Grüneisen Equation of State
The Mie-Grüneisen EOS can be obtained from Equations (45)-(46) with the Hugoniot curve (Equation (45)) and the cold curve (Equation (46)) as the reference curves. The Hugoniot curve is obtained via MD simulation. Based on previous analysis, the Grüneisen coefficient γ in this study is calculated by Equations (32) and (36). The cold pressure and cold energy are calculated from the Born-Mayer potential. Figure 9 shows the Mie-Grüneisen EOS of iron in 3D contour form of SC iron based on the cold pressure and cold energy, with the Hugoniot relations as the reference.
The Mie-Grüneisen EOS is described in the pressure-specific volume-internal energy (P-V-E) space. In the front view, the whole Mie-Grüneisen EOS space presents a concave surface characteristic, as shown in Figure 9. The maximum pressure of EOS calculated from the cold curve or the Hugoniot curve as the reference will be different. The maximum pressure when the Hugoniot relations are utilized as the reference is ~200 GPa higher than that when the cold pressure and cold energy are the references. The results indicate that the reference data have a certain impact on the surface characteristic and values of the Mie-Grüneisen EOS.  The Mie-Grüneisen EOS is described in the pressure-specific volume-internal energy (P-V-E) space. In the front view, the whole Mie-Grüneisen EOS space presents a concave surface characteristic, as shown in Figure 9. The maximum pressure of EOS calculated from the cold curve or the Hugoniot curve as the reference will be different. The maximum pressure when the Hugoniot relations are utilized as the reference is~200 GPa higher than that when the cold pressure and cold energy are the references. The results indicate that the reference data have a certain impact on the surface characteristic and values of the Mie-Grüneisen EOS.

Conclusions
In this study, a micromechanical approach for determining the Mie-Grüneisen EOS parameters of iron under the Hugoniot states was investigated using the MSST method in MD simulation, and the shock Hugoniot, cold curve, melting temperature, and Grüneisen coefficient are also discussed. The simulation results show that the difference between the linear relationships of the shock Huguenot of SC and NC iron was insignificant, and both agree well with the experiment data. Compared with the results at GS 5 nm and 25 nm, the GS effect of iron was not significant at high pressures. The difference between the Born-Mayer potential and the Morse potential in calculating the cold pressure was almost the same as that in calculating the cold energy, but the cold pressure and cold energy in SC iron were a little higher than those in NC iron. As the pressure increased, the Grüneisen coefficient gradually decreased, while the melting temperature gradually increased; the two trends were opposite. The Grüneisen coefficient γ is calculated by Equations (16c) and (18c), which have a better performance in this study. The Mie-Grüneisen EOS presented a concave surface characteristic in 3D contour space; the pressure with the Hugoniot curve as the reference was higher than that with the cold curve as the reference, which indicates that different reference curves have a certain impact on the surface characteristic and values of Mie-Grüneisen EOS. The research results in this paper prove the feasibility of the MSST method and provide a reference for studying the EOS of other metals under high pressures. Materials in practical applications are very complex and contain various impurities and defects, and both impurities and defects affect their properties. We examined pure SC iron and perfect NC iron in this study with the MSST method. In future work, we will consider impurities and defects to study the dynamic response characteristics of materials under more complex conditions. Author Contributions: Conceptualization, investigation and writing-original draft preparation, Y.W.; methodology and software, Y.W. and X.Y.; resources, F.W. and J.D.; supervision, X.Z. and H.C.; project administration and funding acquisition, X.Z. All authors have read and agreed to the published version of the manuscript.
Funding: This research was funded by the foundation of the National Key Defense Laboratory for Computational Physics through the project "Study on the constitutive model of metal materials with damage effect under strong dynamic load", grant number HXo2020-74.

Institutional Review Board Statement: Not applicable.
Informed Consent Statement: Not applicable.

Data Availability Statement:
The data that support the findings of this study are available from the corresponding author upon reasonable request.