Scale effects on hardness during molecular dynamics simulations of nanoindentation
Abstract (After obtaining results)
Several quantities can be evaluated during molecular dynamics simulations of nanoindentation, such as force-displacement curves, and especially dislocation microstructures and densities. Building up on the work by Christoph Begau on the free energy of dislocation microstructures during nanoindentation this thesis will look further into size effects on the evolution of defect structures and densities and the relation to the hardness.
Molecular dynamic simulation of displacement controlled nanoindentation in the (111)-free surface of a Nickel sample will be performed. By properly defining the size of the plastic zone, i.e. restricting it to a small volume under the indenter, the evolution of the Geometrically necessary dislocation (GND) and Statistically stored dislocation (SSD) densities will be analysed for different ranges of indentation depths to a maximum depth of 0.6r, where r is the indenter radius. In addition, the relationship between the hardness and the indentation depth will be investigated to probe the assumption that it is constant.
To this effect the slope of the hardness-vs-displacement curve will be evaluated after relaxations will be taken at certain holding points and the dislocation structures at these points will be analysed.
Talk about MD simulations, Nickel, our study and nanoindentation.
Hardness is a materials resistance to localized plastic deformation. In order to measure the mechanical properties such as hardness or young's modulus of a material, various techniques have been developed in which nanoindentation is one of them and has played a great role in the past years even until date. Nanoindentation simply involves applying a prescribed load to an indenter in contact with a material or specimen. Several parameters can be obtained from this such as the indentation depth, indirect measurement of the contact area - that is, the area of contact between the indenter and the specimen, which is one distinguishing feature of most nanoindentation tests, stress on the material or under the indenter tip just to name a few. Different forms of indenter tips have been applied during nanoindentation such as spherical, Rockwell, Conical, Berkovich, Knoop and Vickers. Each indenter tip has its own effect on the material it is applied on also each material response to an indenter tip is different (Anthony C. Fischer-Cripps, Nanoindentation, 2002, Springer). It is stated that the popularity for spherical indenters are increasing, due to the fact that this type of indenter provides a smooth transition from elastic to elastic-plastic contact . For example, from (MTS test sheet for spherical tip) this is an advantage over Berkovich. Anthony C. Fischer stated that this suitable for measuring soft materials and for replicating contact damage in in-service conditions .
Investigating the deformation mechanisms inside a crystal during plastic deformation can be done by molecular dynamic simulations of nanoindentation. The atomic arrangement or crystalline structure of a material is important to determine its behaviour . Molecular dynamic simulation of nanoindentation enables us to view the evolution of dislocations and the microstructure of this dislocations. In this study, the method of Molecular dynamics is applied for better understanding of the mechanisms that occur during nanoindentation, the scale effects on hardness i.e. the effect of indenter radius, system size, temperature and the volume of the plastic zone under the indenter tip and the evolution of GNDs and SSDs during nanoindentation using Nickel.
Scale effects: clarify them- size of the box, radius of indenter
Motivation for the present study
Nanoindentation experiments and simulation studies have been carried out with Aluminium single crystal to understand the dislocation network and evolution of dislocation densities during this simulation using IMD code, performing this work with Nickel is to broaden our knowledge and see the effects of nano-indentation and the response of different materials surfaces during nano-indentation.
Large-scale MD simulations of nanoindentation in Nickel are performed in order to gain a better understanding concerning the ongoing mechanism of dislocation generation during the initial phases of plastic deformation and the following evolution. Therefore, the objectives of this work are:
To define the plastic zone (right volume and scale relationship) under the indenter
To analyse the evolution of GND and SSD densities during this simulation
To see the effect of temperature on the evolution of the dislocation densities and structures.
To determine the hardness as a function of indentation depth
Organization of thesis
This chapter (Chapter 1) gives a brief introduction of the present work with its primary objectives and a brief look at the MD simulations of the nanoindentation performed using Nickel. Chapter 2 gives a literature review on MD simulations background, various works on Nanoindentation with the effects of Spherical indenter, piling-up and sinking-in effects seen in some materials, brief look at dislocation generation during nanoindentation and, GNDs and SSDs and the plastic zone. Chapter 3 is about MD modelling, features of ITAP IMD and the simulation code and all. Chapter 4 results and discussion. Chapter 5 conclusion of this study and future work.
Molecular Dynamic Simulation
Computer simulations are said to be carried out in the hope of understanding the properties of assemblies of molecules in terms of their structure and the microscopic interactions between them . Michael Allen distinguished two main families of simulation technique as molecular dynamics (MD) and Monte Carlo (MC) where MD simulations has an obvious advantage over MC in the fact that it gives a route to dynamical properties of the system: transport coefficients, time-dependent responses to perturbations, rheological properties and spectra. MD Simulations model the motion of some group of particles by solving the classical equations of motion. While MC methods are much simpler than MD as they do not require the calculation of molecular forces, they do not yield significantly better statistics than MD in each amount of computing time .
Nanoindentation is the most frequently used technique for measuring thin film properties such as Young's modulus and hardness . Understanding and getting information about the mechanical hardness of materials at the nanoscale can be achieved by nanoindentation. Considering crystalline materials, it is stated that nanoindentation can be used to study defect nucleation and propagation events, which are detected by discontinuities in the load-depth relationship . The interest in the experimental and theoretical parts of the events that occur has grown i.e. which are said by Ma, Li, et al., to be sudden displacement-burst excursions during load-controlled nanoindentation of relatively dislocation-free metals.
Fig: Nanoindentation with a spherical indenter. Image taken from (27)
The first pop-in event which are said to have been observed by Gane and Bowde in 1968 is often identified as the initiation of dislocation nucleation, and thus the transition from purely elastic to elastic/plastic deformation . Indentation techniques measure the mechanical properties of thin films and other materials e.g., elastic constants, yield stress, and hardness by measuring the force or load on the indenter tip as a function of the tip displacement during indentation and retraction .
Nanoindentation can also be used to obtain the size effects in metals which was studied by George Z. Voyiadjis in  and stated that these size effects are governed by the behavior of dislocations and their interactions with other types of defects such as grain boundaries. In this case, the plastic zone size divided by the contact radius is not a constant factor and varies as the indentation depth changes and the hardness decreases as the indentation depth increases for conical indenters .
Spherical Indenter and Size Effects
Weichao Guo et al., stated that the mostly used indenters are of two types, the indenters with a revolution surface which includes spherical and conical indenters, and the pyramid indenters which include three-sided Berkovich and four-sided Vickers . Spherical indenters are commonly used for nanoindentation because they provide a smooth transition from elastic to elastic-plastic contact. Experiments to see indentation size effect was studied with different spherical indenter tips were used to show that the indentation hardness H increases with an increase of ratio a / R, (where a is the contact radius and R is the indenter radius) which means larger penetration depth showed larger indentation hardness .
At the same time, H increases with a decrease of the indenter radius R. Numerical Simulations using a spherical indenter have shown to have close results with literature data as described in . Hongwei Zhao et al., investigated the influences of crystal orientations and temperature on the hardness of silicon using a spherical diamond indenter. Results showed that the hardness of monocrystalline silicon becomes smaller with temperature influence, but adhesiveness increases as temperature increases . This was also a result for the Young's modulus and hardness from Te-Hua Fang et al., after studying the effects of temperature on the atomic-scale nanoindentation process . Considering thin-films, atomic simulation of nanoindentation of Ni thin Film were carried out in  to study the effects of indenter velocity and radii on hardness. Both were seen to have negligible effects on the hardness value. Dupont et al., also performed spherical indentation on Ni thin wires to observe crystal plasticity and concluded that the hardness of the nanowires decreases as the sample diameter decreases (25).
READ ref 14 again for more info on Nickel
Piling-Up and Sinking-In
During indentation into an elastic material, the surface of the material is typically drawn inwards and downwards underneath the indenter and sinking-in occurs . When the contact involves plastic deformation, the material may either sink in, or pile up around the indenter. Pile up patterns were studied by (Wang et al., 2004) using a conical indenter for nanoindentation where the pile up patterns on the surface showed four-, two-, and six-fold symmetry, respectively and depend on the plastic work-hardening of the indented material and temperature . Conclusions were made that the piling-up frequently occurs in soft materials whereas the sinking-in occurs in hard materials   and studying the influence of the phenomenon of piling up onto the values of hardness and Young's modulus by nanoindentation tests with spherical/conical tips , results showed that the higher the penetration depth, the larger the relative pile-up effect. Molecular dynamic simulations of nanoindentation of Nickel substrates were carried out in different crystallographic directions to prove that the crystalline nickel materials demonstrate the pile-up phenomenon at maximum indentation depth from nanoindentation . This work should be able to compare the sink-in behavior of Nickel during nanoindentation.
Dislocation Generation During Nanoindentation
Using the Nye tensor analysis, the complete information of the dislocation network becomes available in a representation that matches the conceptual idea of a dislocation network and offers the possibility for further processing . The analysis method to extract dislocation networks introduced as followed is in principle independent of the Nye-tensor data, as already published in Begau et al. .
Dislocation generation is also something to investigate during Nanoindentation. The investigation of dislocation generation and interaction using single-crystal copper and the EAM potential was carried out by C.Begau et al., where interest was placed in understanding how initial dislocations multiply during the early stages of plastic deformation . Also, the use of a skeletonization algorithm to analyze the multiplication of dislocations. Obtaining the dislocation networks include the resultant burgers vectors per atoms which enable the creation of dislocation networks where the skeletonization method also plays a role . This master thesis report will also try to observe the dislocation networks and evolutions using the RBVs during the post processing which enhance skelotonization during nanoindentation of Nickel at different temperatures and indenter radius.
Geometrically Necessary Dislocation Densities and SSDs
The dislocation densities and the geometrically necessary dislocation densities can be computed. Begau stated that the total dislocation density ''_total in an arbitrarily shaped volume v is the total length of all dislocations lines within v (22) as shown in equation
In descriptions of dislocations in crystals, dislocations could be separated into two different categories, as described by  as geometrically-necessary dislocations arising from the geometry of the indenter as described by Anthony Fischer in (1) and by M.F Ashby as simple geometric relationship exists between the deformation and the number of dislocations needed to accommodate it with minimum internal stress due to geometrical constraints of the crystal lattice, and statistically-stored dislocations (SSDs) which evolved from random trapping processes during plastic deformation or arising for statistical reasons. The geometrically necessary dislocations take the form of circular dislocation loops as shown in figure 2;
Fig: Geometrically necessary dislocations in the plastic zone created by a conical indenter. Image taken from (1)
Yu Gao et al., stated 'the region filled with dislocations with increasing indentation depth during nanoindentation is known as the plastic zone' . It was referred to a hemispherical region with a radius equal to the contact radius of the indenter .
The study of the area below the indenter is also a topic to investigate during nanoindentation. Yu Gao et al., studied the size of the plastic zone generated during the nanoindentation in fcc (Cu and Al) and bcc (Fe and Ta) metals by a spherical indenter , the extent of the plastic zone was seen to be substantially larger before the retraction of the indenter. With this assumption, the plastic zone under the indenter during the Nanoindentation of Nickel will be investigated and evaluated.
Procedure for Molecular Dynamic Simulations
Taking a system of N particles, generating trajectories of this system can be done by molecular dynamic simulations by integration of the Newtons equations of motion with specifying the interatomic potential and boundary conditions. MD can be classified into classical MD and ab initio but the focus here is on the former. Taking a system with N-particles and a potential energy E, the atomic trajectories can be found b solving the Newton equations of motion shown in equation;
Where m is the mass and
Molecular Dynamic Simulation Potentials
Simulations involving larger systems require the use of empirical potentials. Assumptions are made and approximations regarding the nature of interaction between atoms and the cost for the evaluation of the energy and forces is reduced so that simulations with millions or billions of atoms are possible. There are different classes of potentials but the focus here is on the Embedded atom potentials which is described in more details in (subsection).
The embedded atom method is currently a common technique used in molecular-dynamics computer simulation of metallic systems . It is said to provide a good description of the interatomic forces in the system, particularly for fcc metals. The interatomic potentials in metals and models forces between atoms can be calculated. First principle DFT calculations led Daw and Baskes (ref) to propose their formula:
Where is the repulsive pair potential and is the embedding energy of atom I with electron density where it is expressed as
General steps of Molecular Dynamic Simulations
a) Initial position and velocity is assigned to every atom.
b) Force on all the atoms are computed with interatomic potential known.
c) Newton's equation of motion is integrated to obtain trajectory of atoms. They are repeated
again, and again until desired property of all the atoms of system is obtained.
Fig: Molecular dynamics Simulations flowchart (26)
Periodic Boundary conditions
The system size that can be simulated with MD is very small compared to real molecular systems . Unwanted artificial boundaries are seen to be present in a system of particles and to avoid this, periodic boundary conditions are introduced which can introduce artificial spatial correlations in too small systems. The use of periodic boundary conditions in our simulations gives us advantages in the sense that if the actual particle should happen to move out of the simulation cell, the image particle in the image cell opposite to the exit side will move in and become the actual particle in the simulation cell. This prevents the loss or creation of particles.
Fig: Periodic boundary conditions showing the original simulation cell colored yellow. The image was taken from 
To keep computations to a manageable level, a force cutoff distance beyond which particle pairs simply do not see each other can also be introduced. The cutoff distance should be set to be less than half of the simulation cell dimension to avoid a particle interacting with its own image.
Some reasons why MD is important are;
Complete control over input, initial and boundary conditions
Detailed atomic trajectories
Several hundred particles are sufficient to simulate bulk matter
Unified study of all physical properties i.e. One can study chemical properties and reactions which are more difficult and will require using quantum MD, or an empirical potential
The major challenge in MD simulations is in the time scale, because most of the processes and experimental observations are at or longer than the time scale of a milli-second.
Itap Molecular Dynamic (IMD)
Developed by the late J''rg Stadler at the Institut f''r Theoretische und Angewandte Physik (ITAP), IMD is a software for classical molecular dynamics simulations which is seen to support several types of interactions such as pair potentials, EAM potentials for metals, Stillinger-Weber and Tersoff potentials for covalent systems, and Gay-Berne potentials for liquid crystals(ref). It is possible to perform simulations of various thermodynamic ensembles also to perform deformation in the process .
As stated in , several simulation options are supported by IMD which allow the following;
shear and deform a sample during the simulation
apply extra forces on certain atoms
fix or rigidly move certain atoms
constrain the mobility of certain atoms
simulate crack propagation
simulate shock waves
simulate laser ablation
compute correlation functions
calculate free energy (http://imd.itap.physik.uni-stuttgart.de/features.html)
IMD runs on the following types of hardware architectures :
Single processor workstations
Workstation clusters with MPI parallelization
Massively parallel supercomputers with MPI parallelization
Multiprocessor SMP machines with OpenMP parallelization
Supercomputers and clusters with several SMP nodes, with MPI parallelization between nodes, and OpenMP parallelization within a node. (http://imd.itap.physik.uni-stuttgart.de/features.html).
The following thermodynamic ensembles are available in IMD as retrieved from (http://imd.itap.physik.uni-stuttgart.de/features.html):
Microcanonical NVE ensemble
Canonical NVT ensemble using a Nos''-Hoover thermostat
Canonical NVT ensemble using a simple Anderson thermostat
Canonical NPT ensemble using a Nos''-Hoover thermostat for both the temperature and the pressure control; there are two variants, one for isotropic volume rescaling, and one which can scale the three axes independently.
Several MD integrators for the relaxation of a structure
Integrators and Relaxators
The integrator is the routine which moves the atoms, depending on the current forces and velocities. The integrators and relaxator used in this work are discussed briefly below. Retrieved from (http://imd.itap.physik.uni-stuttgart.de/userguide/ensembles.html).
The NVE integrator is the constant volume, constant energy ensemble. As stated 'this is the most reliable integrator and thus the main work horse' . No particular parameters are needed, except perhaps the initial temperature start temp, if not all atoms have their velocities given in the configuration file. Retrieved from (http://imd.itap.physik.uni-stuttgart.de/userguide/ensembles.html).
Constant volume, constant temperature ensemble, with Nos''-Hoover thermostat. This is an ensemble with external temperature control where the temperature can actually vary during the simulation (starttemp and endtemp). A further required parameter is the time constant of the temperature control (thermostat), tau_eta. Retrieved from (http://imd.itap.physik.uni-stuttgart.de/userguide/ensembles.html).
If the scalar product of the global force and momentum vectors (containing the force and momentum components of all atoms) is negative (momentum goes "uphill" in the potential landscape), all momenta are reset to zero. The momenta are also reset to zero if the atoms are getting too fast, i.e., the kinetic energy per atom exceeds the value of the parameter glok_ekin_threshold. The glok relaxator is used during the first steps of the simulation where the Nickel structure was relaxed (http://imd.itap.physik.uni-stuttgart.de/userguide/ensembles.html).
Simulation studies Using IMD
For nanoindentation in this present study, the IMD command used to create the structure is:
java -jar JmakeConfig.jar create.param
Indent file which contains all information of the atoms
imd_mpi_eam_glok_nbl_mono_extpot_nye_ada -p indent.param
imd_mpi_eam_nvt_nbl_mono_extpot_nye_ada -p indent.param
Performing the simulations
First the work of Christoph was taken, cell box of 100nm each was used. Indenter radius of 16nm and indentation depth of 0.5r.
Dislocation analysis in FCC and BCC crystals
Calculation of lattice rotation patterns
Calculation of dislocation densities
Import and rendering of any scalar values per atom
Results and discussion
1. Anthony C. Fischer-Cripps, Nanoindentation, 2002, Springer
2. Master thesis report
3. Ma, Li, et al. 'Effect of the Spherical Indenter Tip Assumption on the Initial Plastic Yield Stress.'Nanoindentation in Materials Science, 2012, doi:10.5772/48106.
4.Kelchner, Cynthia L., et al. 'Dislocation Nucleation and Defect Structure during Surface Indentation.' Physical Review B, vol. 58, no. 17, 1998, pp. 11085'11088., doi:10.1103/physrevb.58.11085.
5. Voyiadjis, George Z., and Mohammadreza Yaghoobi. 'Large Scale Atomistic Simulation of Size Effects during Nanoindentation: Dislocation Length and Hardness.' Materials Science and Engineering: A, vol. 634, 2015, pp. 20'31., doi:10.1016/j.msea.2015.03.024.
6. 'On the Influence of Indenter Tip Geometry on the Identification of Material Parameters in Indentation Testing.' Thesis / Dissertation ETD, 2010. Weichao Guo
7. Suresh, S., Nieh, T.G. and Choi, B.W., (1999). Nano-indentation of copper thin films on silicon substrates. Scripta Materialia, Vol. 41, pp. 951-957.
8. Qu, S., Huang, Y., Pharr, G.M. and Hwang, K.C., (2006). The indentation size effect in the spherical indentation of iridium: A study via the conventional theory of mechanism-based strain gradient plasticity. International Journal of Plasticity, Vol. 22, pp. 1265-1286.
9. Wang, Y., Raabe, D., Kluber, C. and Roters, F., (2004). Orientation dependence of nanoindentation pile-up patterns and of nanoindentation microtextures in copper single crystals. Acta Materialia, Vol. 52, pp. 2229-2238.
10. Lee, H., Lee, J.H. and Pharr, G.M., (2005). A numerical approach to spherical indentation techniques for material property evaluation. Journal of the Mechanics and Physics of Solids, Vol. 53, pp. 2037-2069.
11. Hernot, X., Bartier, O., Bekouche, Y., El Abdi, R. and Mauvoisin, G., (2006). Influence of penetration depth and mechanical properties on contact radius determination for spherical indentation. International Journal of Solids and Structures, Vol. 43, pp. 4136-4153
12. Hongwei Zhao, Peng Zhang, Chengli Shi, Chuang Liu, Lei Han, Hongbing Cheng, and Luquan Ren (2014) 'Molecular Dynamics Simulation of the Crystal Orientation and Temperature Influences in the Hardness on Monocrystalline Silicon.' Journal of Nanomaterials, vol. 2014, 2014, pp. 1'8., doi:10.1155/2014/365642.
13. Fang, Te-Hua, et al. 'Molecular Dynamics Analysis of Temperature Effects on Nanoindentation Measurement.' Materials Science and Engineering: A, vol. 357, no. 1-2, 2003, pp. 7'12., doi:10.1016/s0921-5093(03)00219-3.
14. Nair, Arun K., et al. 'Nanoindentation of Thin Films: Simulations and Experiments.' Journal of Materials Research, vol. 24, no. 03, 2009, pp. 1135'1141., doi:10.1557/jmr.2009.0136.
15. Ju, S.-P., et al. 'The Nanoindentation Responses of Nickel Surfaces with Different Crystal Orientations.' Molecular Simulation, vol. 33, no. 11, 2007, pp. 905'917., doi:10.1080/08927020701392954.
16. Maneiro, M.a. Garrido, and J. Rodriguez. 'Pile-up Effect on Nanoindentation Tests with Spherical-Conical Tips.' Scripta Materialia, vol. 52, no. 7, 2005, pp. 593'598., doi:10.1016/j.scriptamat.2004.11.029.
17. Introduction to Molecular Dynamics Simulation Michael P. Allen
18. M. O. Steinhauser and S. Hiermaier, "A Review of Computational Methods in Material Science: Examples from Shock-Wave and Polymer Physics" International of Journal of Molecular Science (2009), 10, 12, 5135-5216.
19. Imd.itap.physik.uni-stuttgart.de. (2018).''IMD - The ITAP Molecular Dynamics Program. [online] Available at: http://imd.itap.physik.uni-stuttgart.de/features.html [Accessed 22 Feb. 2018].
20. Begau, C., et al. 'Atomistic Processes of Dislocation Generation and Plastic Deformation during Nanoindentation.' Acta Materialia, vol. 59, no. 3, 2011, pp. 934'942., doi:10.1016/j.actamat.2010.10.016.
21. Begau, C., et al. 'A Novel Approach to Study Dislocation Density Tensors and Lattice Rotation Patterns in Atomistic Simulations.' Journal of the Mechanics and Physics of Solids, vol. 60, no. 4, 2012, pp. 711'722., doi:10.1016/j.jmps.2011.12.005.
22. Begau, Christoph. 'Characterization of Crystal Defects during Molecular Dynamics Simulations of Mechanical Deformation.' Bochum, 2012.
23. Ashby, M. F. 'The Deformation of Plastically Non-Homogeneous Materials.' Philosophical Magazine, vol. 21, no. 170, 1970, pp. 399'424., doi:10.1080/14786437008238426
24. Gao, Yu, et al. 'Comparative Simulation Study of the Structure of the Plastic Zone Produced by Nanoindentation.' Journal of the Mechanics and Physics of Solids, vol. 75, 2015, pp. 58'75., doi:10.1016/j.jmps.2014.11.005.
25. V. Dupont and F. Sansoz Molecular dynamics study of crystal plasticity during nanoindentation in Ni nanowires
26. Masumeh Foroutan, S. Mahmood Fatemi Molecular Dynamics Simulations of Systems Containing Nanostructures
27. Zhao, Hongwei & Zhang, Peng & Shi, Chengli & Liu, Chuang & Han, Lei & Cheng, Hongbing & Ren, Luquan. (2014). Molecular Dynamics Simulation of the Crystal Orientation and Temperature Influences in the Hardness on Monocrystalline Silicon. Journal of Nanomaterials. 2014. 1-8. 10.1155/2014/365642.
Zhao, Hongwei, et al. 'Molecular Dynamics Simulation of the Crystal Orientation and Temperature Influences in the Hardness on Monocrystalline Silicon.'Journal of Nanomaterials, vol. 2014, 2014, pp. 1'8., doi:10.1155/2014/365642.
...(download the rest of the essay above)