Molecular Dynamics Simulation of Tensile Behavior on Ceramic Particles Reinforced Aluminum Matrix Nanocomposites
Heng Gu1, *, Jiao Jiao Wang2, Zhongwei Li3
1Rexa, Inc. West Bridgewater, MA, USA
2School of Materials Science and Engineering, East China University of Science and Technology, Shanghai, China
3Department of Naval Architecture and Marine Engineering, University of New Orleans, New Orleans, LA, USA
To cite this article:
Heng Gu, Jiao Jiao Wang, Zhongwei Li. Molecular Dynamics Simulation of Tensile Behavior on Ceramic Particles Reinforced Aluminum Matrix Nanocomposites. International Journal of Materials Science and Applications. Vol. 5, No. 3, 2016, pp. 151-159. doi: 10.11648/j.ijmsa.20160503.16
Received: June 10, 2016; Accepted: July 4, 2016; Published: July 7, 2016
Abstract: The mechanical properties and interfacial structures for aluminum matrix composites reinforced by nanometer-sized SiC-β particles has been studied using molecular dynamics (MD) simulation. The modified embedded atom methods, was implemented to describe the atomic interactions. The molecular model undergoes an annealing MD simulation from 300 K to 1000 K to reach its minimum energy point. Tensile tests were performed with periodic boundary conditions. The stress-strain relationship has been studied and elastic constants were predicted as well. The results were compared with those given by continuum-based finite element analysis (FEA) together with the experimental data available in the literatures. It showed that both the elastic modulus and yield stress were further strengthened due to the presence of the nano-particles. Also, it was found that the existing SiC nano-particles have an effect on the initial arrangement of Al atoms in such a manner: Al atoms were inclined to aggregate around the particle surface. Aluminum concentrations were also observed inside the SiC particles close to the surface. The depth of hybridization is uniform and planar.
Keywords: Molecular Dynamics, Aluminum Matrix Composites, Silicon Carbide Nano-particles
In the past decades, metal matrix composites (MMCs) have emerged as a class of materials capable of advanced structural, aerospace, automotive, electronic, thermal management, and wear applications [1-2,3]. The selective matrix alloys, such as Aluminium [4,5], have been predominantly employed for the production of MMCs reinforced by nanoparticles due to the fact that it is superior in terms of improved mechanical and thermal properties (including specific strength and modulus, elevated temperature stability, higher thermal conductivity, and controlled coefficient of thermal expansion) [6,7,8].
The large surface-volume ratio of the nano-sized inclusions, significant improvements in the material properties have been observed at the nanoscale and has become active in recent years [9,10,11,12]. Such mechanical property improvements have resulted in major interest in nanocomposite materials in numerous automotive and general industrial applications . These include potential for utilization as mirror housings on various vehicle types, door handles, engine covers and intake manifolds and timing belt covers. More general applications currently being considered include usage as impellers and blades for vacuum cleaners, aircraft parts, power tool housings, mower hoods and covers for portable electronics such as mobile phones and pagers [14,15,16].
Although substantial technical challenges remain in fabrication, several methods have been developed, e.g. Yang and Li  presented an inexpensive fabrication method for achieving uniform distribution of nano-sized SiC particles in molten aluminum alloy A356, employing the ultrasonic nonlinear effects, known as transit cavitation and acoustic streaming . A total of 20% hardness increasing was measured with a 2.0 vol% nanoparticle addition in aluminum alloy. Gierlotka et al.  used high-pressure infiltration methods to synthesize metal-ceramic nanocomposites: this process is done in a toroid high pressure (7.7 GPa) and high temperature (up to 20000°C). The nanoporous matrix is prepared by compacting nanopowders of high-hardness materials such as Al2O3, SiC, or diamond et al. [20,21,22]. The infiltrating material can be any substance with the melting temperature (at a given pressure) below the processing temperature. When the infiltrant melts the pressure forces it to fill the pores in the matrix. The resulting composite contains a continuous network of solidified injected material with embedded grains of the ceramic powder [23,24].
Due to the small length scale, neither classical continuum mechanics nor traditional numerical methods, such as FEA, are sufficient in predicting the material properties [25,26,27]. The atomistic interaction between a ceramics and a metal is of great significance in the improvement of the bulk material properties of the composites. Traditional numerical methods do not take this into account. Molecular dynamics (MD) could well describe the atomistic interaction between two condensed-phase materials and provides a proper numerical approach to characterization of nanocomposites . However, the high computational cost associated with large-scale MD simulations limits its use beyond the nanoscale whereas the use of micromechanical homogenization techniques in a multiscale framework can lead to the propagation of large errors across length scales . Therefore, efforts have been made to maximize computational efficiency and maintain high accuracy while investigating the material behavior at the micro-, meso and macro-scales [30,31,32]. A significant amount of research has been conducted in this field to develop a multiscale model to characterize material properties [33,34]. Techniques that rely on discrete material distribution and simplified analyses have been developed. Underhill and Doyle  developed a spring-bead model in which each bead represents a large segment of molecular chains, to model the mechanical response of a polymer. Furthermore, Zhang et al. [36,37] developed a spring-bead based cross linked network model which can simulate larger polymer chains and microstructures. In Zhang’s work, MD simulation was implemented to provide necessary information to construct the model. The network model was found to be more accurate than micromechanical approximation techniques for characterizing amorphous materials. This result indicates the MD simulation is important in development of multiscale model of composite materials.
In this work, a molecular model is developed, which is characterized by a representative volume element (RVE) of Al matrix with a silicon carbide (SiC) particle embedded inside. The interface between the Al matrix and SiC particle was built up in a proper way to avoid the defects produced by the misfit of the lattice constants in the initial structure. The nanocomposite model was built at a relatively high temperature (1000 K) at which the aluminum is in the liquid state, which was cooled down to room temperature (300 K). The calculated results have proved that there is an apparent enhancement in material properties including tensile strength and toughness. The morphological features were also examined: the aluminum atoms in the neighborhood of the particle align with the shape of the inclusion. Al atoms are easy to be captured to cluster around carbon-rich sites, so hybridization could be observed beneath the surface of SiC within a certain depth.
2. Basic Theories in Molecular Dynamics
Calculations in molecular dynamics (MD) are based on the energy minimization. An appropriate interatomic potential plays a key role in MD simulations. The potential law between atoms in this study obeys Lennard-Jones potential as follows .
Although this is an approximate potential, it can character the interactions between closed aligned atoms, which is a strong repulsive short range interaction, namely a long range van der Waals attraction. The mathematical form shown in Eq.1 is simple, but it is actually more realistic in its treatment of long range interactions, than the Morse function. However, it only deals with a single pair of atoms. For further-step accurate simulations, many body interactions will be applied. This requires a through understanding of the influence of the local electronic density. In the recent few decades, lots of progress has been made in identifying multi-atom potentials. The most popular theories are: the embedded atom model (EAM), the finnis-Sinclair potentials and the tight-binding potential. Based on different approaches, these models are shown to be equivalent from the mathematical point of view. In this work, the modified embedded atom method (MEAM)  will be employed to describe the intermetallic and metal-ceramic interactions as well as the bonding in the SiC-β particle. It has been shown that MEAM matches well with empirical and ab-initio calculations [38,40]. In the MEAM the total energy E of the system is as the sum of defined energy atoms, with each energy contributions from an embedding function F connecting to the electron density and a pair potential Φ:
Here A is a constant. The pair potential, ф, only depends on the distance between atoms i and j, Rij=. Ec is the cohesive energy. is a reference electron density. The actual angular dependence of the potential is embedded in the definition of . which represents the electron density at atomic site i. The electron density is defined as:
where t is an adjustable parameter which weights the directional dependence, is a spherically symmetric term, (h = 1,2,3) is the angular dependent terms of the electron density as
is the electron density from a neighboring atom at a distance R from the atom of interest. The α, β and γ variables are summed over the three directional components of the distance vector between atoms i and j. 
3. Model of Molecular Dynamics Simulation
As shown in Fig. 2 (a), pure aluminum is characterized by face-centered cubic (FCC) periodic crystal structure. Silicon carbide (Fig. 2 (b) has the same FCC crystal structure, but in the compound, each silicon atom is bonded to four carbon atoms. The covalent bonds between silicon and carbon are very strong. Thus the bulk silicon carbide presents high stiffness. As shown in Fig. 3, the model consists of a cube of aluminum atoms with a center-placed spherical inclusion of silicon carbide (red part). The edge length of the cube is 4.050nm and the diameter of the sphere equals 1.305nm.
The molecular model containing 4,719 atoms was set up and schematically shown in Fig. 4. The aluminum was aligned in a matrix with a block SiC-β embedded inside to ensure both materials have common crystalline orientation. Each inner surface of the aluminum block is parallel to the corresponding outer surface of the sub-block. The interface geometries were created by assigning a proper value to the interfacial gap. This interface separation value t was chosen to be 1.575Å, at which the lattice mismatch and the influence of residual strain in the initial model so that the whole structure turns out to be periodic.
The Si and C atoms were placed within a spherical region whose geometric center coincides with that of the blocks. The radius of the sphere, r, was set to be 6.525Å, which is equal to half of the edge length of the sub-block. The redundant Si and C atoms outside the sphere (the grey part) were removed. Afterwards, a second spherical surface, sharing the same geometric center and with a radius of 8.10Å (r+t), was created. The vacuum space outside the surface was filled in with aluminum atomsl the atoms are displayed with their van der Waals radius, and each the nanoparticle can be approximated to be a sphere, with the diameter of 1.305nm.
4. Energy Minimization
Molecular mechanics consider the molecule to be a collection of atoms held together by simple elastic or harmonic forces. The MD simulation uses the parameters of force fields, a mathematical expressions employing a combination of internal coordinates and bond terms (bond distances, bond angles, torsions, etc.), to describe the dependence of the energy of a molecule on the coordinates of the atoms in the molecule, and non-bond terms to describe the van der Waals and electrostatic interactions between atoms.
For condense-phase materials such as ceramics and metals considered here, the force field given by the Condensed-phase Optimized MolecularPotential (COMPASS) was used, whose development focuses on high accuracy.
The initial temperature of the simulation box was set to be 0 K using the theory of molecular mechanics. After this calculation the system was equilibrated. Since the minimization methods can only optimize the atomistic model to the closest local minimum point (Fig. 7) , so this state of equilibrium largely depends on the starting structure, from where the system energy falls to a local minimum point.
The most straightforward way to find the global minimum is to do a systematic conformational search. However, it is time-consuming and infeasible for systems with large number of molecules. As a result, simulated annealing is commonly used to find the global minimum. In the simulated annealing procedure, the simulation commences at a high temperature where most attempted (moves) are accepted and the temperature is then slowly reduced. The transition probability is reduced in concert so that the average energy of the sampled points in the space also diminishes. At the conclusion of the annealing, the resulting point will, in general, be near the global energy minimum for the system.
The atoms were initially assigned with velocities following a Maxwell-Boltzmann distribution around the target temperature (300 K). Then the interatomic potential (MEAM), a dynamics simulation was performed in our model for 0.001 ps, with the isothermal canonical ensemble NVT (constant volume-constant temperature) using Nosė-Hoover method as thermostat to control the temperature. The calculation started at an initial temperature of 300 K which would be raised up to the mid-cycle temperature 1100 K and then cooled down to 300K again. The temperature increment was 50 K. There were 500 time steps in maximum for integration. The model underwent 15 cycles.
5. Results and discussion
5.1. Interfacial Morphology
Interfacial structure and chemistry are important factors influencing the properties of material. The final configuration obtained in this work is shown in Fig. 9 (a). In the far field of the matrix, the atom arrangement is still in a good agreement with crystal structure of pure aluminum. It is clearly shown that the morphology of atoms in the center of the SiC-β nanosphere do not change too much due to the strong bonds between silicon and carbon atoms, but in its neighboring region to the aluminum matrix original lattice constants are severely distorted. The atoms are clustered and follow the structure of the particle. This actually increases the effective size of the nanoparticle. In order to study the influence of size effect, the same MD calculation was performed on a larger system with particle of 2.18nm in diameter (Fig. 9(b)). This result is similar to that in Fig. 9(a), with clearer interfacial layer observed.
The chemical composition of the Al/SiC interfacial structure using two different annealing conditions has been studied by Arsennalt R. J. et al . At a higher temperature above 700°C, Al4C3 formed rapidly to form a planer and clean interfacial layer. The formation rate of Al4C3 is quite fast so that the diffusion of Al into SiC is prevented. Bermudez V.  investigate the the physical and electronic structure of the interface between Al and SiC, employing spectroscopic methods. His results revealed that the aluminum atoms form islands randomly distributed over the surface of SiC particle and Al aggregates at Carbon-rich points and reacts with Carbon to form Al4C3. Theoretical quantum mechanics studies were also carried out: Arsenault R. J., Li S. and Jena P.  used the intermediate neglect of differential overlap (INDO) method to study adhesion at SiC/Al interfaces and proved that the bond strength between SiC and Al could be 1.5 times stronger than the bond between Al and Al. Rao and Jena  performed a self-consistent quantum mechanical calculation to study the Al-Al, Al-Si, Al-C interactions using unrestricted Hartree-Fock theory. Their results show that the binding energy of Al-C is significantly higher than that of Al-Al and Al-Si. When an interface is formed between Al and SiC, Al will be bonded with the exposed C atoms, which accounts for the annealing study in .
As shown in Fig. 6, after minimizing the energy at 0 K, the SiC phase looks chemically stable. No hybridization has been observed, so we did not obtain any evidence that diffusions or penetrations of Al atoms into SiC compound took place. However, in Fig. 9, the Al atoms became clustered around the sphere after annealing above 1000K, forming a layer of Al islands distributed randomly over the surface of the SiC particle. The Al layer is non-uniform but turns out to be patching on the SiC surface. It is noticed that Al atoms are more likely to deposit around the exposed C atoms. Former quantum chemical calculations well explained it: the Al-C interaction is much stronger than Al-Al, Al-Si, so Al atoms are more possible to be captured by C. The Al concentration of about 50% of the Al matrix is seen to be present 2-3Å deep into the SiC. The hybridization band of Al and SiC is quite neat and of uniform thickness, since when it goes beyond 3Å from the interface into the particle, the Al concentration falls sharply to very small magnitude close to 0. This matches the experimental result that the formation of a planar Al4C3 layer prevents the diffusion of Al into SiC as well as the dissolving of SiC . As shown in Fig. 10, Si and C concentrations diminish very quickly when it goes out of the interface.
5.2. Elastic Constants
Assigning boundary conditions to the periodic model, a simulated tension test was carried out. By fixing the x=0 surface, a strain increment of 0.1% is applied to another surface which is perpendicular to x-axis at 300 K. Then the system is relaxed by minimizing the energy. There are various theories for the calculation of mechanical properties. The bulk modulus of the metals and alloys was calculated by using the method proposed in . The elastic constants were calculated at the experimental lattice constants from the difference in the total energies of the distorted and undistorted lattices. In Hu’s work, the elastic constants of Gold Nanorods Produced are determined by Seed Mediated Growth . The methodology is to use a single-point energy calculation to obtain the second derivatives of the lattice energy with respect to the lattice parameters and the atomic coordinates. The energy expression:
where Uo is the equilibrium energy and are the strain. When the system is at an energy minimum (that is, all first-order derivatives of the lattice energy are zero), the second-order derivative term can be used to calculate the components Cijkl of the stiffness tensor
A finite element analysis was also performed using the classical continuum theory by assuming the composites simply consist of two phases (the interfacial layer would not be taken into account). The effective volume fraction of the SiC particle was taken to be identical to that of the molecular RVE. The FE analysis used the same boundary conditions and dimensions as those for the MD model. Both Al and Al6061 are used for the ductile matrix: E = 68.0GPa, ν = 0.33 for pure Al; E = 68.3GPa, ν = 0.33 for Al6061; E = 401GPa, ν = 0.22 For silicon carbide. The stress-strain relations obtained from MD simulation and FE analyses are compared (Fig. 11 and Fig. 12).
The results indicate that the Young’s moduli estimated based on MD simulations are 8-10% higher than the results obtained from FEA. The yield strength from MD is 15% higher than that from FEA. This MD results match well with the reported results in which a 2.0vol% nanoparticle addition in aluminum alloy results in a 20% higher stiffness . Hardening occurs when the strain is higher than 4.0%. The stress is raised up to 290MPa. The simplified FE analysis in which there existed only two phases does not take into account of the molecular structure of the metal/ceramic interface. There is high stress concentration at the interface between Al and SiC. The response of the composites to the applied loading is also a function of the spatial distribution of the inclusions, given that the volume fraction does not change. The elastic constants from molecular dynamics present a small anisotropy in the orthogonal directions. This is mainly due to the geometrical imperfection of the spherical SiC nanoparticle induced by the approximate discretization: the particle surface is not perfectly spherical and its radius is only a nominal value. FE analysis cannot predict the elastic constants accurately because the continuum model regards the material as uniform, while the density changes along the radial distance from the interface due to the rearrangement of the atoms in the matrix and the mass conservation. The mechanical reinforcement of the composites revealed in the current MD study can be explained by the strong chemical bonds between aluminum and carbon atoms, as illustrated in previous experiments and quantum mechanics calculations . Surely, the calculation was based on the assumption that there is no defect in the material and the inclusion is perfectly bonded to the matrix in reality.
In this work, molecular dynamics simulations have been performed on an Al-SiC nano-composite system using the interatomic potential. By employing annealing technique, the system reach an equilibrium state at the room temperature. Based on the principle of energy minimization, the elastic constants for the composites were estimated and the SiC/Al interfacial structure was also investigated. The results reveal that the strong bonding energy between aluminum and carbon can significantly strengthen the composites and an effective interface can be obtained. The results match well with existing experiments and FE simulations.
The thoughtful discussion with Dr. J. Zhang and his advices are gratefully appreciated. The authors also wish to thank Dr. K. Li for his contributions and Dr. X. Li for providing experimental data. The computations are supported by supercomputer facility and Laboratory of Molecular Simulation (LMS) at Texas A&M University, which are thankfully acknowledged.