Al solid–liquid interfaces: A molecular dynamics study

Al solid–liquid interfaces: A molecular dynamics study

Computational Materials Science 174 (2020) 109489 Contents lists available at ScienceDirect Computational Materials Science journal homepage: www.el...

3MB Sizes 0 Downloads 17 Views

Computational Materials Science 174 (2020) 109489

Contents lists available at ScienceDirect

Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci

The orientation dependence of liquid ordering at α-Al2O3/Al solid–liquid interfaces: A molecular dynamics study ⁎

R. Yana, W.Z. Suna, S.D. Maa, T. Jinga, , H.B. Dongb,

T



a Key Laboratory for Advanced Materials Processing Technology, Ministry of Education, School of Materials Science and Engineering, Tsinghua University, Beijing 100084, China b School of Engineering, University of Leicester, Leicester LE1 7RH, UK

A R T I C LE I N FO

A B S T R A C T

Keywords: Molecular dynamics simulation ReaxFF Solid–liquid interface Liquid ordering Temperature effect

In this study, the evolution of ordering in liquid Al adjacent to α-Al2O3 substrates with (0 0 0 1) and (1120 ) surface orientations was examined using a reactive force field (ReaxFF) and molecular dynamics (MD) simu-





lation. It is found that atoms from the liquid would first form a layer at the (0001) or (1120 ) interfaces to fulfill the crystal periodicity. Atoms in this layer adopt the in-plane structure of the substrate and lack mobility. Once the crystal periodicity is fulfilled, liquid layers near the (0 0 0 1) interface exhibit an in-plane structure similar to −

the (1 1 1) plane in face-centered cubic (FCC) Al, while liquid layers near the (1120 ) interface show a hexagonal structure. Distances between neighboring liquid layers do not show a dependence on temperature, but the extent of ordering in liquid Al increases with decreasing temperature.

1. Introduction Numerical and experimental studies have proven that a transitional region exists near the interface between a substrate and the adjacent liquid [1–6]. There are two forms of ordering within the transitional region. One is out-of-plane layering, which describes damped density oscillations perpendicular to the solid–liquid interface. The other is inplane ordering, which characterizes the arrangement of atoms in the ordered liquid layers. Since a deeper understanding in the structure of solid–liquid interface would promote the development of a wide range of fundamental science and technological processes, such as wetting [7,8], nucleation [9–12], crystal growth [4,13,14] and additive manufacturing [15], the study of ordering at solid–liquid interfaces has attracted a lot of attention in recent years. As for out-of-plane layering, Broughton et al. [16] simulated the FCC (1 0 0), (1 1 0) and (1 1 1) solid–liquid interfaces using a modified Lennard-Jones potential. They found that the distance between ordered layers in liquid adopts the interlayer spacing d of the substrate, and when d spacing is close to the nearest neighbor distance of bulk liquid, out-of-plane layering persists in a longer distance. Their finding was in agreement with results obtained in later simulations [2,17–20] using more realistic potentials. With the aid of high resolution transmission electron microscopy (HRTEM) and a numerical phase retrieval method, Kauffmann et al. [21] made a quantitative analysis of the α-



Al2O3(0001)/Al solid–liquid interface, showing that the distance between the first ordered Al layers is quite similar to the d spacing of the α-Al2O3(0001) planes, which agrees with Broughton’s finding. Later, Gandman et al. [22] used the aberration corrected transmission electron microscopy to investigate the α-Al2O3/Al solid–liquid interfaces. Unlike Kauffmann’s result, they found that the distances between the ordered Al layers near the α-Al2O3(0001) surface vary between the d spacing of the α-Al2O3(0001) planes and the nearest neighbor distance in bulk liquid Al. Moreover, Gandman et al. found that for the α− Al2O3(1012 )/Al solid–liquid interface, the distance between the ordered layers mimics the nearest neighbor distance of the bulk liquid instead of the d spacing of the substrate. This phenomenon was also observed at the Cu/Pb and Al/Pb solid–liquid interfaces [23,24]. In addition to these two cases, significant interfacial alloying was observed at the Fe (1 1 1)/Li solid–liquid interface [25], resulting in to no apparent out-ofplane layering at 800 K and 1100 K. The measure of in-plane ordering is much more difficult than that of out-of-plane layering either in terms of experiments or numerical simulations. As a result, much fewer quantitative analyses of in-plane ordering have been conducted. Most studies [2,16,19,20,22,24–26] have shown that the in-plane structure of liquid layers is similar to the counterpart of the adjacent substrate. However, Hashibon et al. [27] simulated the interfaces between liquid Al and BCC substrates and found that the original BCC (1 0 0) surface of the substrate transformed

Corresponding authors. E-mail addresses: [email protected] (T. Jing), [email protected] (H.B. Dong).

https://doi.org/10.1016/j.commatsci.2019.109489 Received 11 October 2019; Received in revised form 6 December 2019; Accepted 15 December 2019 0927-0256/ © 2019 Elsevier B.V. All rights reserved.

Computational Materials Science 174 (2020) 109489

R. Yan, et al.

this force field to investigate the effect of carbon coating on the oxidation of Al nanoparticles. Their simulation results are in agreement with those from experiments. Hong et al. also reported that the melting point of aluminum nanoparticle is ~800 K, but they did not give the precise data. Our preliminary study found that the melting point of bulk Al is 860 ± 5 K. This value will help to determine the temperature range used in this study. However, when the α-Al2O3 sample was put together with the liquid Al sample, an interdiffusion phenomenon occurred at the interface: O atoms in the α-Al2O3 diffused into the liquid and Al atoms from the liquid diffused into the α-Al2O3. This phenomenon has been reported by Streitz et al. [36]who examined the solid–solid interface between Al (1 1 1) and α-Al2O3 (0001) using an ES + potential and by Zhang et al. [37] who investigated both the solid–solid and solid–liquid interfaces between Al and α-Al2O3 (0001) using the ReaxFF potential. Since experiments [5,21,22] and AIMD simulations [30,31] all showed that the α-Al2O3 would not undergo substantial structure change, this interdiffusion phenomenon is certainly a flaw of the ReaxFF. Considering the significance of MD simulations of α-Al2O3/Al solid–liquid interfaces, we made a compromise to avoid the interdiffusion phenomenon. In our simulations, atoms in the α-Al2O3 substrates were fixed and only atoms from the liquid Al were allowed to move when the substrates and the liquid were put together. The method of pinning atoms in the substrate was employed in previous studies of ordering at solid–liquid interfaces [17,18] and impressing results have been obtained.

to an FCC (1 0 0) surface due to adsorption of liquid Al atoms at the interface. The experiments conducted by Gandman et al. [22] also revealed that the in-plane periodicity in liquid Al adjacent to α− Al2O3(1012 ) plane differs from that in the substrate. Moreover, prefreezing of Pb layers has been reported at the Cu(1 1 1)/Pb solid–liquid interface [23], with a lattice spacing 33% larger than that of the Cu substrate and a rotation angle of 6° relative to the underlying Cu(1 1 1) substrate. It can be seen from the above that the ordering phenomena at solid–liquid interface are very complicated and a correlation between the liquid ordering and the substrate structure is still lacking. Furthermore, due to the difficulties encountered in experiments and simulations, very few temperature points have been used, resulting in the evolution of the structure at solid–liquid interfaces with temperature less addressed. In this study, we examine the solid–liquid interfaces between liquid Al and − α-Al2O3 (0001) and (1120 ) planes at different temperatures using MD simulation and the latest version of ReaxFF for Al-O description − [28,29]. MD simulations of α-Al2O3(0001)/Al and α-Al2O3 (1120 )/Al solid–liquid interfaces are important for three reasons. First, so far experimental studies [5,11,21,22] about solid–liquid interfaces have focused on the α-Al2O3/Al system and most of them [5,11,21] have only examined the α-Al2O3(0001)/Al interface, investigations of other interfaces are needed. Moreover, although these experiments have successfully extracted quantitative information about the out-of-plane layering, they have not given the structure within the layers. To our knowledge, no MD simulations about ordering at α-Al2O3/Al solid–liquid interfaces have yet been reported. Therefore, quantitative information about the out-of-plane layering at the α-Al2O3(0001)/Al interface obtained in our MD simulation could be compared with results from experiments and ab initio molecular dynamics (AIMD) simulations to evaluate the performance of the potential. If good matches are achieved, then the in-plane-ordering at the (0001) interface and the − results for the (1120 ) interface predicted by our MD simulation are more convincing and can serve as important complements to previous experimental studies. Second, just as mentioned above, the influence of liquid ordering on subsequent nucleation process is worth exploring. Although AIMD simulations [30–32] have been used to study the αAl2O3(0001)/Al solid–liquid interface, they are not efficient enough to handle this problem due to the requirements for large simulation size and long simulation time. Since MD simulations make it possible to simulate bigger systems and to run longer time, our MD simulation of − the (0001) and (1120 ) interfaces could be further extended to examine the effect of liquid ordering on solidification.

2.2. Interface construction All MD simulations were carried out using the LAMMPS program [38] with a time step of 1 fs. A Nosé–Hoover thermostat and a Nosé–Hoover barostat were used to control temperature and pressure, respectively, with relaxation times being 100 fs and 1000 fs. Visualization was performed by the OVITO software [39]. In this study, two substrate surface orientations, i.e., α-Al2O3(0001) − and α-Al2O3(1120 ) were examined. The simulation began with creating bulk α-Al2O3 samples at 1000 K. For the (0001) orientation, the z axis was kept along the α-Al2O3 [0 0 0 1] direction, while the × and y axes − − − along the [1210 ] and [1010 ] directions, respectively. The dimensions along the × , y and z axes were 8a, 4 3 a, 2c, respectively, where a = 4.76 Å, c = 13.0 Å [40] are the lattice constants of α-Al2O3. The initial configuration of the created α-Al2O3 sample with the (0001) orientation is shown in Fig. 1(a) (top view) and Fig. 1(b) (side view). In an α-Al2O3 unit cell, O atoms are hexagonal close packed with the stacking sequence of ABAB, and the arrangements of O atoms in layer A and layer B are the same; while Al atoms locate at octahedral sites formed by the O atoms. The stacking sequence of Al atoms is ABCABC, where the A or B or C layer is marked with a red rectangle in Fig. 1(b) and its top view is given in Fig. 1(c). In fact, the A (or B or C) layer consists of two Al layers. These two Al layers have the same in-plane structure, which is shown in Fig. 1(d), but they don’t align with each other vertically and thus form the pattern shown in Fig. 1(c) when seen from the z axis. − To create the bulk α-Al2O3 sample with the (1120 ) surface orienta− − tion, the × , y and z axes were along the α-Al2O3 [1100 ], [0001], [1120 ] directions, respectively, and the corresponding dimensions were 4 3 a, 3c, 6a. The initial configuration of the created α-Al2O3 sample with the − (1120 ) surface orientation is given in Fig. 2(a) (top view) and Fig. 2(b) (side view). It can be seen from Fig. 2(b) that along the z axis only a single Al layer, which is marked out by a black rectangle, occurs between two adjacent O layers. The arrangement of Al atoms in this single layer is shown in Fig. 2(c). After getting the initial configurations, periodic conditions were applied in all directions. The α-Al2O3 samples were then equilibrated in NPT simulations for 105 steps with the target temperature and target pressure being 1000 K and 0 bar, respectively. During the NPT runs,

2. Simulation details 2.1. Potential selection Choosing a suitable potential for the modelled system is vital for the success of MD simulations. However, it is very difficult to develop a good potential for simulating α-Al2O3/Al solid–liquid interfaces because this system contains both metallic and ionic bonds and because charge transfer effects at the interface must be taken into account. In the process of potential selection, we considered three potentials which involve dynamic charge equilibration: the CTIP potential proposed by Zhou et al. [33], the COMB3 potential developed by Choudhary et al. [34] and the ReaxFF by Duin et al. [28,29]. These potentials were first used to model an α-Al2O3 crystal in the MD simulation. We found that the α-Al2O3 crystal modelled by the CTIP potential or by the COMB3 potential could not keep stable above 500 K, while the crystal modelled by the ReaxFF behaved very well even above 1000 K. Therefore, the ReaxFF was chosen for our simulation. The ReaxFF incorporates bond order concept and dynamic charge equilibration [35], which are essential for modelling complex phenomena occurring at the interfaces between liquid Al and α-Al2O3 substrates. Hong et al. [29] has utilized 2

Computational Materials Science 174 (2020) 109489

R. Yan, et al.

Fig. 1. Atom arrangements in the created αAl2O3 sample with the (0001) surface orientation. Red and blue balls represent Al and O atoms, respectively. (a) Top view of the whole sample. (b) Side view of the whole sample, and the black and red rectangles denote regions containing single and double Al layers, respectively. (c) Top view of the double Al layers. (d) Top view of the single Al layer.

solid–liquid interface the termination of the α-Al2O3 surface could be Al termination or O termination, namely, a layer of Al or O atoms could be exposed on the substrate surface and come into contact with the liquid. − In this study, for both of the (0001) and (1120 ) only the O termination was considered. And vacuum regions were added below the α-Al2O3 samples and above the Al samples to avoid possible effects of pressure on the solid–liquid interfaces. periodic conditions were only applied to × and y directions in following simulations. Fig. 3 shows the coexistence system with the α-Al2O3(0001)/Al solid–liquid interface.

potential energies, temperatures, pressures and system sizes were recorded. After 104 steps, all these quantities began to slightly fluctuate around their own equilibrium values. For each configuration, a snapshot with sizes that matched the average sizes was chosen to construct the corresponding solid–liquid coexistence system. To obtain liquid Al samples, an Al crystal containing 4000 atoms was first melt in a NPT run lasting for 105 steps, with the target temperature and target pressure being 2000 K and 0 bar, after which the temperature was switched to 1000 K and 105 more steps of NPT run was carried out. The obtained configuration was then used in two different runs. For each of the two runs, the cross section of the obtained Al sample was first adjusted to that of the α-Al2O3 substrate. Then the final structures were put in two NPzAT runs, which fixed the cross section of the Al samples and allowed the Al samples to change their sizes in the z axis. All of the NPzAT runs lasted for 2 × 105 steps to make the Al samples reach equilibrium. Lastly, the equilibrated α-Al2O3 samples were combined with their corresponding Al samples respectively, with an initial gap of 1 Å between the solid part and liquid part. It should be noted that at the

2.3. Interface characterization Ordering in liquid at temperatures around its melting point is of great importance. As the melting point of Al reproduced by the ReaxFF is 860 ± 5 K, three temperatures were utilized in this study: 1000 K, 900 K and 800 K. After constructing the interfaces, the α-Al2O3 substrates in all coexistence systems were fixed and only Al atoms were allowed to move. The coexistence systems were first equilibrated at 1000 K for 105 steps using NVT simulations, followed by 2 × 105 steps



Fig. 2. Atom arrangements in the created α-Al2O3 sample with the (1120 ) surface orientation. Red and blue balls represent Al and O atoms, respectively. (a) Top view of the whole sample. (b) Side view of the whole sample, and the black rectangle denotes a single Al layer. (c) Top view of the single Al layer. 3

Computational Materials Science 174 (2020) 109489

R. Yan, et al.

where i refers to the atom type (Al in α-Al2O3, O in α-Al2O3 or Al from the liquid), Nzi is the number of i atoms locating between z-Δz/2 and z+Δz/2 and the angled brackets indicate the time average of the 200 recorded configurations, and L x and L y are dimensions of the system along the x and y axes, respectively. Layers in the transitional region near the interface can be defined as regions between minima of the density profile. The in-plane ordering is examined by 2D density profile and 2D structure factor map. Each layer was divided into small areas with dimensions of 0.25 Å×0.25 Å in the xy plane. The distribution of atoms in the xy plane was used to obtain the instantaneous 2D density profile:

N→ r ρxy (→ r)= Δx Δy

(2)

where N→ r indicates the number of all atoms locating in the area with − positon r in a recorded configuration, Δx and Δy are the dimensions of the small areas along the x and y axes, respectively, which are both 0.25 Å The instantaneous 2D density profile was then used to obtain the 2D structure factor map by averaging the 200 instantaneous 2D density profiles, which is defined as −



Fxy (k ) = 〈 |ρxy (k )|2 〉

(3)



where ρxy (k ) is the Fourier transform of the instantaneous 2D density profile. The dynamical properties near the solid–liquid interface were investigated by calculating diffusion constants in the layers. For each of the layers, atoms locating within the layer centered at z at time t0 were assigned to the same group, then the average mean-square displacement (MSD) per atom 〈 [r j (t ) − r j (t0 )] 〉z for this group was computed every 20 steps for 500 times. Therefore, the MSD is function of time t. Then the diffusion constant in this layer was derived by determining the limit of

Fig. 3. Initial configuration of the simulation system with the α-Al2O3(0001)/ Al solid–liquid interface. Red, blue and green balls represent Al atoms in αAl2O3, O atoms in α-Al2O3 and Al atoms in liquid Al, respectively.

D (z ) =

1 d 〈 [r j (t ) − r j (t0 )] 〉z 6 dt

(4)

To improve the statistics, the MSD in each layer was calculated from 10 independent time origins separated by 104 steps. And to make a comparison with diffusion constant in bulk liquid, for each of the interface diffusion constant in the liquid region far from the interface was also calculated. 3. Results and discussion 3.1. Out-of-plane layering −

Fig. 5 shows the density profiles for the (0001) and (1120) interfaces at 1000 K, 900 K and 800 K, with z = 0 aligning to the center of the outmost fixed O layer (denoted by 0 in Fig. 5) of the substrate. It can be seen clearly that whatever the substrate orientation is, several liquid Al layers form near the interface and these layers decay rapidly with distance from the interface. Fig. 5 also shows that the density profiles depend on temperature, namely, the lower the temperature is, the higher the peaks in the liquid are. Apart from these similarities, each interface also has its own unique characteristics. For the (0001) interface, some Al atoms from the liquid penetrate into the substrate and form a layer (denoted by −1 in Fig. 5(a)) below the outmost O layers. This can be explained that in the α-Al2O3 substrate only two thirds of the octahedral sites formed by neighboring O atoms are occupied by Al atoms and thus Al atoms from the liquid easily occupy the left octahedral sites due to the interdiffusion tendency. Therefore, the −1 layer will not be further examined in the rest of this study. On the other side of the outmost O layer, the first Al layer splits into two sublayers (denoted by S1 and S2 in Fig. 5(a)), which is different from the other layers formed in the liquid. This finding is in agreement

Fig. 4. Temperature profile of the NVT runs for interface characterization.

of sampling. The systems were then cooled from 1000 K to 900 K using 105 steps, and 3 × 105 more steps were performed at 900 K, with the first 105 steps for equilibrium and the last 2 × 105 steps for sampling. The cooling and sampling processes were repeated until sampling at 800 K was carried out. Fig. 4 illustrates the change of temperature during the whole simulations. In the sampling stage, the systems were divided into uniform bins with thickness (Δz ) of 0.08 Å along the z direction. Configurations were recorded every 1000 steps for 200 times. To characterize the out-ofplane layering, the distribution of atoms along the z direction in every recorded configuration was used to calculate the density profile, which is defined as

ρi (z ) =

〈Nzi 〉 L x L y Δz

(1) 4

Computational Materials Science 174 (2020) 109489

R. Yan, et al.



Fig. 5. Density profiles for the (a) α-Al2O3(0001)/Al and (b) α-Al2O3(1120 )/Al solid–liquid interfaces, with red dot lines representing best fittings of density profiles

(

at 800 K to an exponential function: ρ (z ) = Aexp −

z ̂ I¾

) + B.

with studies of Kang et al. [30] and Fang et al. [32] who both examined the same interface using AIMD simulations. Later we will see the inplane structure of these two sublayers also differentiates from that in − the following layers. At the (1120 ) interface more detectable peaks can be found with decreasing temperature. Some atoms from the liquid move to the substrate surface and form a layer (denoted by 1) with its peak height being much higher than those of layers immersed in the liquid. After the first layer, a gap appears, which totally separates the first and second layer. This is not observed at the (0001) interface and will be discussed in the section below. 3.2. In-plane ordering More information can be drawn from these layers’ in-plane structures. Fig. 6 shows the 2D density profiles and 2D structure factor maps for different layers near the (0001) interface at 800 K. The first layer’s 2D density map is similar to that of the double Al layers in the substrate (Fig. 1(c)), however, it is also obvious that some spots in the 2D density map are brighter than the others. Meanwhile, its structure factor map exhibits both regular spots corresponding to long range ordering and rings representing bulk liquid. As the first layer splits into two sublayers, the in-plane structures of these two sublayers were examined separately. It becomes clear that the brighter spots in the first layer’s 2D density profile mainly come from the S1 sublayer while the darker spots mainly from the S2 sublayer. Moreover, while the S1 sublayer is quite ordered as regular spots are still clear even far from the center of its structure factor plot, some vacancies can be seen in the S2 sublayer’s 2D density map and spots only appear in the center region of its structure factor plot. This means that the S2 sublayer is more mobile than the S1 sublayer. An AIMD calculation carried out by Kang et al. [30] also shows that these two sublayers adopt the structure of the double Al layers in the substrate (Fig. 1(c)) and vacancies exist exclusively in the S2 sublayer. The agreement between our result and the work of Kang et al. shows that the ReaxFF is capable of simulating the α-Al2O3/Al solid–liquid interfaces. When it comes to the layer 2, a hexagon appears in the center of its structure factor plot, combined with several rings. Three basic repeated cells are marked out in its 2D density profile using red lines and red dots. Obviously, layer 2 does not adopt the in-plane structure of the substrate, instead, its in-plane structure is similar to the (1 1 1) plane of FCC Al. This is quite impressing as experiments [12,41] have found that Al(1 1 1) // α-Al2O3(0001) is the preferred orientation relationship at the α-Al2O3/Al solid–solid interface. In addition, Ma et al. [31] performed an AIMD calculation, in which they used an orientation order parameter to analyze the atom arrangement in the liquid layers. They found that the liquid Al layer adjacent to αAl2O3(0001) is weakly six-fold symmetrical, agreeing with the result

Fig. 6. 2D density profiles (left column) and structure factor profiles (right column) for selected layers at 800 K near the α-Al2O3(0001)/Al solid–liquid interface.

5

Computational Materials Science 174 (2020) 109489

R. Yan, et al.

second ring resembling a hexagon. These two features are totally absent in the structure factor map of layer 1, indicating that the in-plane structure of layer 2 is different from that of layer 1. The in-plane structure here is also different from that at the (0001) interfaces, as no atoms could be found in the center of the hexagon. The 2D density profile of layer 3 becomes obscurer and no more features except rings could be found in its structure factor map. This indicates that in-plane ordering does not exist in layer 3. Note that for clarity all the 2D density profiles and 2D structure factor maps in this section are from simulation runs at 800 K. In-plane structures in these layers at 1000 K and 900 K are similar to those at 800 K, but corresponding features in the 2D density profiles and structure factor maps are less obvious. The effect of temperature on the liquid ordering will be evaluated in section 3.4. Now we use the outmost O layer in the substrate as a reference and examine atom arrangements on the liquid side. Interestingly, for either − the (0001) or (1120 ) interface, nearest to the outmost O layer, atoms from the liquid will first form an layer which adopts the in-plane structure of the substrate. In other words, atoms near the interface will fulfill the periodicity of the substrate first. Once the periodicity of the substrate is fulfilled and can no longer continue (no O atoms in the liquid in this study), subsequent liquid layers will have their own inplane structures as a result of the competition between the substrate structure and the bulk liquid structure. 3.3. Interfacial diffusion profiles The dynamical property in the interfacial region is characterized by the diffusion profile, which is shown in Fig. 8. The diffusion profiles show an obvious dependence on temperature as expected. At each temperature, the diffusion coefficients in the bulk liquid calculated from the two systems are nearly identical, demonstrating that the parameters we used for this calculation are appropriate. Moreover, diffusion coefficients in the bulk liquid at 1000 K agree pretty well with experimental value [42], which further validates the performance of the ReaxFF. At the (0001) interface the diffusion constant in the S1 sublayer is nearly zero, smaller than that in the S2 sublayer. This difference is in accordance with the result of their 2D density maps and structure factor plots. However, the diffusion constant in the S2 sublayer is still much smaller than those in the subsequent layers. This means that the motion of atoms in the S2 sublayer is reduced compared with those in the layers followed. From the second layer, the diffusion coefficient begins to increase monotonically and almost reaches the − value in bulk liquid in the fifth layer. For the (1120 ) interface, the diffusion coefficient in the first layer is essentially zero, meaning that although this layer is formed by atoms from the liquid, it is actually a solid layer. Atoms in the second layer are more mobile, resulting in a

Fig. 7. 2D density profiles (left column) and structure factor profiles (right −

column) for selected layers at 800 K near the α-Al2O3(1120 )/Al solid–liquid interface.

given here. The in-plane ordering fades quickly and almost disappears in layer 3 because no spots except rings could be seen in its structure factor profile. − It can be found from Fig. 7(a) that layer 1 at the (1120 ) interface − mimics the structure of the single Al layer in the α-Al2O3(1120 ) substrate (Fig. 2(c)). The regular spots in the structure factor plot indicate layer 1 is quite stable and the in-plane ordering is strong. Interestingly, layer 2 exhibits a hexagonal structure similar to that of graphene. Three basic repeated cells are marked out with red hexagons in its 2D density map. The hexagonal structure in layer 2 could be demonstrated by two features in its structure factor plot: one is the six spots forming a hexagon in the first (or innermost) ring and the other is the distorted



Fig. 8. Diffusion profiles near the (a) α-Al2O3(0001)/Al and (b) α-Al2O3(1120 )/Al solid–liquid interfaces. For clarity, the corresponding density profiles were also attached. Error bars represent the standard deviations evaluated from the 10 origins. 6

Computational Materials Science 174 (2020) 109489

R. Yan, et al.



Fig. 9. Distances between neighboring layers for the (a) α-Al2O3(0001)/Al and (b) α-Al2O3(1120 )/Al solid–liquid interfaces at 1000 K, 900 K and 800 K. Distance at layer i is calculated by using z i + 1 − z i , where z i is the center of layer i. Error bars represent the standard deviations evaluated from the 5 density profiles and each of which was computed from 40 consecutive recorded frames.

performed from the second layer, while two ways of fitting for the (0001) interface was utilized: one is to fit from the S1 sublayer, and the other is to fit from layer 2. Superscripts in table 1 tell from which layer fittings began. Our calculation shows the extent of ordering at the interface increases with decreasing temperature regardless of substrate orientation. It is also found that at 1000 K the calculated correlation lengths for the − (0001) and (1120 ) interfaces are 3.43 A ̂ ± 0.25 Å and 3.10 A ̂ ± 0.17 Å, respectively, agreeing very well with the experimental values at a slightly higher temperature 1023 K. It should also be noted that for the (0001) interface correlation lengths extracted from the second layer are much greater than those extracted from the S1 sublayer. This is not surprising because although fitting could be done from the S1 sublayer, the in-plane structure of the is different from that of subsequent layers, namely, the ordering extent decreases when including the. The significant difference between values fitted from the S1 sublayer and from the layer 2 highlights the importance of the choice of layers when deriving the correlation length from a density profile. This may be the reason for the discrepancies between the calculated correlation lengths at 900 K and the experimental value at 943 K for the (0001) interface as the identification of the first liquid layer is difficult in experiments.

much larger diffusion coefficient in this layer. After a period of monotonic increase, the diffusion coefficient in the sixth layer is almost equal to that in the bulk liquid. 3.4. Temperature effect on liquid ordering The temperature effect is first investigated by comparing distances between neighboring liquid layers at different temperatures. To get the layering distance, the 200 recorded frames were used again and five more density profiles were obtained by averaging every 40 of the recorded frames. The five density profiles were then utilized to extract the centers of the layer peaks, from which the distances between neighboring layers and their errors were evaluated. It should be noted that since the first layer at the (0001) interface splits into two sublayers, the average value of the peaks of the two sublayers is regarded as the center of the first layer at the (0001) interface. It can be seen from Fig. 9 that although peaks shift to the substrate with decreasing temperature (see Fig. 5), distances between neighboring layers do not undergo noticeable changes with temperature regardless of substrate orientation. In addition, at the (0001) interface distances between neighboring layers are very close to each other, fluctuating around 2.45 Å. Although distances − between layer 2, 3, 4 and 5 at the (1120 ) interface are also about 2.45 Å, distance between the first and the second layer is greater than this value. This larger layering distance may be due to the structure difference between these two layers, which results in a poor wetting of the − liquid to the substrate at the (1120 ) interface. However, because the substrate surface is totally covered by the liquid in our simulation, it is impossible to detect the contact angle here. We will test this hypothesis by cooling the system to lower temperature in our future studies. Gandman et al. [22] investigated the α-Al2O3(0001)/Al and α− Al2O3(1120 )/Al solid–liquid interfaces using HRTEM. They found that distances between neighboring Al layers are 2–3 Å and 2.22–2.66 Å, respectively, agreeing very well with results given here. The temperature effect is further examined by fitting the density

(

z

4. Conclusion In this study we conducted the first MD simulation which examined the evolution of ordering phenomenon in liquid Al adjacent to α-Al2O3 − substrates with (0001) and (1120 ) surface orientations using a ReaxFF potential. The two solid–liquid interfaces were characterized through calculations of density profiles, 2D density profiles, structure factor plot and diffusion profiles. Layering distance and correlation length have also been extracted from the calculated density profiles. The calculated in-plane structure of the first layer at the (0001) interface agrees very well with the AIMD simulation [30]. The diffusion coefficient in bulk Al and correlation length obtained in this study were also compared with available experimental values [21,22,42] and good matches were achieved. These agreements validate the ability of the ReaxFF to simulate the α-Al2O3/Al solid–liquid interfaces although atoms in the αAl2O3 substrate were fixed in the simulation. Our analysis shows that − for either the (0001) or (1120 ) interface, atoms from the liquid would first form a layer to fulfill the crystal periodicity. Once the crystal periodicity is fulfilled, an in-plane structure similar to FCC Al (1 1 1) plane appears at the (0001) interface, while a hexagonal structure − presents at the (1120 ) interface. In-plane ordering exists in less than

)

profile to an exponential function ρ (z ) = Aexp − ̂ + B , where A is a I¾ normalization factor, B is a constant representing the density of bulk ̂ is the correlation length which is widely used to liquid, and I ¾ quantitatively describe the ordering extent at the interface. Table 1 shows values of correlation length for different α-Al2O3/Al solid–liquid interfaces from this work and from experiments carried out by researchers. During the fitting process, we found that if the first layer at − the (1120) interface is taken into account, the fitting would fail because −

it is too high. In this work, the fitting for the (1120) interface was 7

Computational Materials Science 174 (2020) 109489

R. Yan, et al.

Table 1 Values of correlation length (Å) for different α-Al2O3/Al solid–liquid interfaces. Note that the superscript (S1 or 2) means that correlation length was obtained by fitting density profile from that layer. For example, (0001) 2 means that correlation length was obtained by fitting density profile from layer 2 in Fig. 5(a). α-Al2O3 surface

MD calculation (this work)

experiments

800 K

900 K

1000 K

943 K

1023 K

3.88 A ̂ ± 0.23 4.15 A ̂ ± 0.21

3.43 A ̂ ± 0.25 3.10 A ̂ ± 0.17





(1120 ) (0001) S1

4.72 A ̂ ± 0.38 4.99 A ̂ ± 0.24





2.69 A ̂ ± 0.35

2.42 A ̂ ± 0.21

2.34 A ̂ ± 0.09





(0001)















2.50 A ̂ ± 0.32 [22] –

3.01 A ̂ ± 0.28 [21] 2.76 A ̂ ± 0.47 [22]

(0001) −

2 2

(1120 )

three layers at all examined interfaces. Distances between neighboring layers at the two interfaces are not influenced by temperature and are very close to each other except the distance between the first and the − second layer at the (1120 ) interface which is greater than the rest. Ordering in the liquid shows a strong dependence on temperature. For any interface, the correlation length increases with decreasing temperature. It is also found that when extracting correlation length by fitting density profile, the choice of the first layer included in the fitting process has a great influence on the obtained value. As experiments about α-Al2O3/Al solid–liquid interfaces have been limited mainly on the (0001) inter− face, our simulation of the (1120 ) interface provides a good complement to experimental studies.

modelled with COMB3 potential, Comput. Mater. Sci. 155 (2018) 136–143. [11] S. Ma, A.J. Brown, R. Yan, R.L. Davidchack, P.B. Howes, C. Nicklin, Q. Zhai, T. Jing, H. Dong, Atomistics of pre-nucleation layering of liquid metals at the interface with poor nucleants, Commun. Chem. 2 (1) (2019) 1. [12] A.J. Brown, H.B. Dong, P.B. Howes, C.L. Nicklin, In situ observation of the orientation relationship at the interface plane between substrate and nucleus using Xray scattering techniques, Scr. Mater. 77 (2014) 60–63. [13] H. Men, A molecular dynamics study of growth anisotropy in Al melt, Modell. Simul. Mater. Sci. Eng. 24 (1) (2015) 015008. [14] Y. Mishin, M. Asta, J. Li, Atomistic modeling of interfaces and their impact on microstructure and properties, Acta Mater. 58 (4) (2010) 1117–1151. [15] W. Sun, R. Yan, Y. Zhang, H. Dong, T. Jing, GPU-accelerated three-dimensional large-scale simulation of dendrite growth for Ti6Al4V alloy based on multi-component phase-field model, Comput. Mater. Sci. 160 (2019) 149–158. [16] J.Q. Broughton, G.H. Gilmer, Molecular dynamics of the crystal–fluid interface. V. Structure and dynamics of crystal–melt systems, J. Chem. Phys. 84 (10) (1986) 5749–5758. [17] A. Hashibon, J. Adler, M.W. Finnis, W.D. Kaplan, Atomistic study of structural correlations at a liquid–solid interface, Comput. Mater. Sci. 24 (4) (2002) 443–452. [18] H. Men, Z. Fan, Atomic ordering in liquid aluminium induced by substrates with misfits, Comput. Mater. Sci. 85 (85) (2014) 1–7. [19] G.Q. Yang, J.F. Li, Q.W. Shi, L.T. Kong, Structural and dynamical properties of heterogeneous solid–liquid Ta–Cu interfaces: a molecular dynamics study, Comput. Mater. Sci. 86 (2014) 64–72. [20] C. Qi, J. Li, B. Xu, L. Kong, S. Zhao, Atomistic characterization of solid-liquid interfaces in the Cu-Ni binary alloy system, Comput. Mater. Sci. 125 (2016) 72–81. [21] Y. Kauffmann, S.H. Oh, C.T. Koch, A. Hashibon, C. Scheu, M. Ruhle, W.D. Kaplan, Quantitative analysis of layering and in-plane structural ordering at an aluminaaluminum solid-liquid interface, Acta Mater. 59 (11) (2011) 4378–4386. [22] G. Maria, K. Yaron, C.T. Koch, W.D. Kaplan, Direct quantification of ordering at a solid-liquid interface using aberration corrected transmission electron microscopy, Phys. Rev. Lett. 110 (8) (2013) 086106. [23] J.P. Palafox-Hernandez, B.B. Laird, M. Asta, Atomistic characterization of the Cu–Pb solid–liquid interface, Acta Mater. 59 (8) (2011) 3137–3144. [24] Y. Yang, D.L. Olmsted, M. Asta, B.B. Laird, Atomistic characterization of the chemically heterogeneous Al–Pb solid–liquid interface, Acta Mater. 60 (12) (2012) 4960–4971. [25] X. Gan, S. Xiao, H. Deng, X. Li, W. Hu, Orientation dependences of the Fe-Li solidliquid interface properties: atomistic simulations, J. Alloy. Compd. 687 (2016) 875–884. [26] R. Yan, S. Ma, T. Jing, H. Dong, The in-plane structure and dynamic property of the homogeneous Al-Al solid-liquid interface, Metals 8 (8) (2018) 602. [27] A. Hashibon, J. Adler, M.W. Finnis, W.D. Kaplan, Ordering at solid-liquid interfaces between dissimilar materials, Interface Sci. 9 (3–4) (2001) 175–181. [28] S. Hong, A.C.T.V. Duin, Molecular dynamics simulations of the oxidation of aluminum nanoparticles using the reaxff reactive force field, J. Phys. Chem. C 119 (31) (2015) 17876–17886. [29] S. Hong, A.C.T.V. Duin, Atomistic-scale analysis of carbon coating and its effect on the oxidation of aluminum nanoparticles by ReaxFF-molecular dynamics simulations, J. Phys. Chem. C 120 (17) (2016) 9464–9474. [30] J. Kang, J. Zhu, C. Curtis, D. Blake, G. Glatzmaier, Y.-H. Kim, S.-H. Wei, Atomically abrupt liquid-oxide interface stabilized by self-regulated interfacial defects: the case of Al/Al2O3 interfaces, Phys. Rev. Lett. 108 (22) (2012). [31] S. Ma, Y. Rui, J. Tao, H. Dong, Substrate-induced liquid layering: a new insight into the heterogeneous nucleation of liquid metals, Metals 8 (7) (2018). [32] C. Fang, Z. Fan, Prenucleation at the liquid-Al/α-Al2O3 and the liquid-Al/MgO interfaces, Comput. Mater. Sci. 171 (2020) 109258. [33] X.W. Zhou, H.N.G. Wadley, J.S. Filhol, M.N. Neurock, Modified charge transferembedded atom method potential for metal metal oxide systems, Phys. Rev. B 69 (3) (2004). [34] K. Choudhary, T. Liang, A. Chernatynskiy, S.R. Phillpot, S.B. Sinnott, Charge optimized many-body (COMB) potential for Al2O3 materials, interfaces, and nanostructures, J. Phys.-Conden. Matter 27 (30) (2015). [35] H.M. Aktulga, J.C. Fogarty, S.A. Pandit, A.Y. Grama, Parallel reactive molecular dynamics: numerical methods and algorithmic techniques, Parallel Comput. 38 (4–5) (2012) 245–259. [36] F.H. Streitz, J.W. Mintmire, Metal/oxide interfaces: an electrostatics-based model, Compos. Interf. 2 (6) (1994) 473–484.

CRediT authorship contribution statement R. Yan: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Resources, Software, Validation, Visualization, Writing - original draft. W.Z. Sun: Data curation, Validation, Writing - review & editing. S.D. Ma: Validation, Writing review & editing. T. Jing: Funding acquisition, Project administration, Supervision, Writing - review & editing. H.B. Dong: Funding acquisition, Project administration, Supervision, Writing - review & editing. Acknowledgements This work was supported by the National Natural Science Foundation of China (51674153, 51320105003) and the National Key Research and Development Programm of China (2017YFB1103700). The calculations were performed on the TianHe-1(A) at National Supercomputer Centre in Tianjin. References [1] R.L. Davidchack, B.B. Laird, Simulation of the binary hard-sphere crystal/melt interface, Phys. Rev. E 54 (6) (1996) R5905. [2] R.L. Davidchack, B.B. Laird, Simulation of the hard-sphere crystal–melt interface, J. Chem. Phys. 108 (22) (1998) 9452–9462. [3] H. Men Z. Fan Molecular dynamic simulation of the atomic structure of aluminum solid-liquid interfaces Materials Research Express 1(2) 025705 (13 2014 pp. )-025705. [4] M. Guerdane, H. Teichler, B. Nestler, Local atomic order in the melt and solid-liquid interface effect on the growth kinetics in a metallic alloy model, Phys. Rev. Lett. 110 (8) (2013) 086105. [5] S.H. Oh, Y. Kauffmann, C. Scheu, W.D. Kaplan, M. Ruhle, Ordered liquid aluminum at the interface with sapphire, Science 310 (5748) (2005) 661–663. [6] T.U. Schulli, R. Daudin, G. Renaud, A. Vaysset, O. Geaymond, A. Pasturel, Substrate-enhanced supercooling in AuSi eutectic droplets, Nature 464 (7292) (2010) 1174–1177. [7] E.B. Webb, G.S. Grest, D.R. Heine, Precursor film controlled wetting of Pb on Cu, Phys. Rev. Lett. 91 (23) (2003) 236102. [8] J. Aguilarsantillan, The influence of surface anisotropy crystalline structure on wetting of sapphire by molten aluminum, Metallurgical and Materials Transactions A-physical Metallurgy and Materials Science 44(5) (2013) 2299–2306. [9] J.P. Palafox-Hernandez, B.B. Laird, Orientation dependence of heterogeneous nucleation at the Cu–Pb solid-liquid interface, J. Chem. Phys. 145 (21) (2016) 211914. [10] R. Yan, W. Sun, S. Ma, R.L. Davidchack, N.D. Pasquale, Q. Zhai, T. Jing, H. Dong, Structural and mechanical properties of homogeneous solid-liquid interface of Al

8

Computational Materials Science 174 (2020) 109489

R. Yan, et al.

and ionic mobility of alpha-Al2O3 at high temperatures, J. Am. Ceram. Soc. 65 (9) (1982) 460–464. [41] M. Hai, Q. Liu, L. Liu, L. Xin, W. She, P. Zhai, Molecular dynamics simulations of the microstructure of the aluminum/alumina interfacial layer, Appl. Surf. Sci. 324 (2015) 538–546. [42] F. Kargl, E. Sondermann, H. Weis, A. Meyer, Impact of convective flow on longcapillary chemical diffusion studies of liquid binary alloys, High Temperat.-High Pressur. 42 (1) (2013) 3–21.

[37] Q. Zhang, T. Cagin, A. van Duin, W.A. Goddard, Y. Qi, L.G. Hector, Adhesion and nonwetting-wetting transition in the Al/alpha-Al2O3 interface, Phys. Rev. B 69 (4) (2004). [38] S. Plimpton, Fast parallel algorithms for short-range molecular dynamics, J. Comput. Phys. 117 (1) (1995) 1–19. [39] A. Stukowski, Visualization and analysis of atomistic simulation data with OVITOthe open visualization tool, Modell. Simul. Mater. Sci. Eng. 18 (1) (2010). [40] P. Aldebert, J.P. Traverse, Neutron-diffraction study of structural characteristics

9