In this research, we have worked on the brittle fracture of graphene nano-ribbon to explore the behavior of crack propagation at different crack angles. We have performed classical Molecular Dynamics simulations using LAMMPS at ten different crack angles between 0 degrees and 45 degrees, in an increment of 5 degrees to observe the parameters that dominate the crack path. The graphene nanoribbon is loaded in the zigzag direction by pulling it in the armchair direction with a pre-existing crack in the center. We have used OVITO for the visualization of the simulation. AIREBO potential is employed in this work because it is extensively used in the fracture of graphene with different loading conditions and temperatures. The crack path is determined for all ten nanoribbons and the nanoribbon with a crack at 25° turned out to be the weakest because of the sharp crack tip and crack shape. The results are validated with the published results and are in accordance with them

Pathak, Vijay Kumar

English

Michigan Technological University

Michigan Technological University Digital Commons @ Michigan Tech Dissertations, Master's Theses and Master's Reports 2021 STUDYING THE EFFECTS OF INITIAL CRACK ANGLE ON THE CRACK PROPAGATION IN GRAPHENE NANO-RIBBON THROUGH MOLECULAR DYNAMICS SIMULATIONS Vijay Kumar Pathak Michigan Technological University, vpathak@mtu.edu Copyright 2021 Vijay Kumar Pathak Recommended Citation Pathak, Vijay Kumar, "STUDYING THE EFFECTS OF INITIAL CRACK ANGLE ON THE CRACK PROPAGATION IN GRAPHENE NANO-RIBBON THROUGH MOLECULAR DYNAMICS SIMULATIONS", Open Access Master's Report, Michigan Technological University, 2021. https://doi.org/10.37099/mtu.dc.etdr/1210 Follow this and additional works at: https://digitalcommons.mtu.edu/etdr  Part of the Applied Mechanics Commons, Computer-Aided Engineering and Design Commons, and the Semiconductor and Optical Materials Commons STUDYING THE EFFECTS OF INITIAL CRACK ANGLE ON THE CRACKPROPAGATION IN GRAPHENE NANO-RIBBON THROUGH MOLECULARDYNAMICS SIMULATIONSByVijay Kumar PathakA REPORTSubmitted in partial fulfillment of the requirements for the degree ofMASTER OF SCIENCEIn Mechanical EngineeringMICHIGAN TECHNOLOGICAL UNIVERSITY2021© 2021 Vijay Kumar PathakThis report has been approved in partial fulfillment of the requirements for the Degreeof MASTER OF SCIENCE in Mechanical Engineering.Department of Mechanical Engineering - Engineering MechanicsReport Advisor: Dr. Susanta GhoshCommittee Member: Dr. Dibakar DattaCommittee Member: Dr. Vijaya V. N. Sriram MalladiDepartment Chair: Dr. William W. PredebonDedicationTo my mother, teachers and friendswho didn’t hesitate to criticize my work at every stage - without which I would neitherbe who I am nor would this work be what it is today.ContentsList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viiAcknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viiiAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ix1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1 Graphene . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.2 Molecular Dynamics . . . . . . . . . . . . . . . . . . . . . . . . . . 22 Theory and Practice . . . . . . . . . . . . . . . . . . . . . . . . . . . 42.1 LAMMPS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52.2 Griffith’s brittle fracture . . . . . . . . . . . . . . . . . . . . . . . . 73 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 83.1 GNR with a pre-crack oriented at 0° . . . . . . . . . . . . . . . . . 93.2 GNR with a pre-crack oriented at 5° . . . . . . . . . . . . . . . . . 103.3 GNR with a pre-crack oriented at 10° . . . . . . . . . . . . . . . . . 11iv3.4 GNR with a pre-crack oriented at 15° . . . . . . . . . . . . . . . . . 123.5 GNR with a pre-crack oriented at 20° . . . . . . . . . . . . . . . . . 133.6 GNR with a pre-crack oriented at 25° . . . . . . . . . . . . . . . . . 143.7 GNR with a pre-crack oriented at 30° . . . . . . . . . . . . . . . . . 153.8 GNR with a pre-crack oriented at 35° . . . . . . . . . . . . . . . . . 163.9 GNR with a pre-crack oriented at 40° . . . . . . . . . . . . . . . . . 173.10 GNR with a pre-crack oriented at 45° . . . . . . . . . . . . . . . . . 183.11 Virial stress computation in LAMMPS . . . . . . . . . . . . . . . . 183.12 Stress-Strain curve . . . . . . . . . . . . . . . . . . . . . . . . . . . 203.13 Validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 213.14 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 223.15 Future work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23vList of Figures1.1 Honeycomb structure of monolayer graphene. Figure taken from [1] 21.2 Schematic of a Molecular Dynamics simulation. Figure taken from [2] 32.1 GNR inside the simulation box with a pre-crack oriented at 0°. The xand y axes correspond to the horizontal and vertical axes respectively. 52.2 GNR with a pre-crack oriented at 45°. . . . . . . . . . . . . . . . . 63.1 Crack propagation in GNR with a pre-crack at 0°. . . . . . . . . . . 93.2 Crack propagation in GNR with a pre-crack at 5°. . . . . . . . . . . 103.3 Crack propagation in GNR with a pre-crack at 10°. . . . . . . . . . 113.4 Crack propagation in GNR with a pre-crack at 15°. . . . . . . . . . 123.5 Crack propagation in GNR with a pre-crack at 20°. . . . . . . . . . 133.6 Crack propagation in GNR with a pre-crack at 25°. . . . . . . . . . 143.7 Crack propagation in GNR with a pre-crack at 30°. . . . . . . . . . 153.8 Crack propagation in GNR with a pre-crack at 35°. . . . . . . . . . 163.9 Crack propagation in GNR with a pre-crack at 40°. . . . . . . . . . 173.10 Crack propagation in GNR with a pre-crack at 45°. . . . . . . . . . 183.11 Stress-strain curve for all GNRs . . . . . . . . . . . . . . . . . . . . 20viList of Tables3.1 Calculated Kc values corresponding to their crack angles . . . . . . 213.2 Calculated σc√a0 values corresponding to their crack angles . . . . 21viiAcknowledgmentsI would like to first thank my advisor, Dr. Susanta Ghosh, for formulating the researchtopic for me and for helping me with everything, by providing necessary resources andconstant guidance. I have learned a lot from Dr. Ghosh and his passion for learninghas certainly ignited an eternal knowledge-seeking flame within me.I would like to thank Dr. Dibakar Datta for teaching me how to perform MD simu-lation and the theory behind it. Your constant support and knowledge-sharing havehelped me a lot in getting this work done.I would also like to thank Dr. Sriram Malladi for his support as a committee mem-ber. I would like to thank my research group members, especially Upendra Yadavfor his initial guidance regarding the use of computational resources and ShashankPathrudkar for his constant feedback on this work as it will continue.In addition, I would like to thank Michigan Technological University for all the re-sources I have utilized to complete my report. I have used the HPC facility (Superior)at Michigan Tech for all my simulations, so I want to thank them separately.I would like to thank all those who have helped me learn, understand and appreciatethis subject as well as those who helped me with LATEX.viiiAbstractIn this research, we have worked on the brittle fracture of graphene nano-ribbonto explore the behavior of crack propagation at different crack angles. We haveperformed classical Molecular Dynamics simulations using LAMMPS at ten differentcrack angles between 0 degrees and 45 degrees, in an increment of 5 degrees to observethe parameters that dominate the crack path. The graphene nanoribbon is loadedin the zigzag direction by pulling it in the armchair direction with a pre-existingcrack in the center. We have used OVITO for the visualization of the simulation.AIREBO potential is employed in this work because it is extensively used in thefracture of graphene with different loading conditions and temperatures. The crackpath is determined for all ten nanoribbons and the nanoribbon with a crack at 25°turned out to be the weakest because of the sharp crack tip and crack shape. Theresults are validated with the published results and are in accordance with them.ixChapter 1Introduction1.1 GrapheneIn this age, when there is a huge requirement for synthesized materials of specificmaterial properties, Graphene has become popular because it is the strongest knownmaterial. Not only it has excellent mechanical properties [3], but it also has very highthermal conductivity [4] (potentially a high-demand material in the semiconductorindustry) and optical properties [5].Graphene is a single layer of carbon atoms arranged in a honeycomb manner. Thishoneycomb shape enables a material to have minimal density and relatively highout-of-plane compression and shear properties.[6]1Figure 1.1: Honeycomb structure of monolayer graphene. Figure takenfrom [1]This is one of the reasons why it has Young’s Modulus of 1 Terapascal and intrinsicstrength of 130 Gigapascals.[3] It is a super-material and can be used in multipleindustries like biomedical, composites, electronics, sensors, energy, etc.1.2 Molecular DynamicsMolecular Dynamics is the study of movements of atoms/molecules which is per-formed through a computer simulation. It is employed for studies in the area ofmaterial science, biophysics, chemical physics, etc. It is the calculation of trajecto-ries of atoms using Newton’s equations of motion. All interactions between differentatoms in a system are given in an interatomic potential, which is used to calculate2the forces between atoms and potential energies of the atoms.[7]Figure 1.2: Schematic of a Molecular Dynamics simulation. Figure takenfrom [2]In the initialization part, first, the atomistic model (geometry) is set up and then themodel is defined through an input file that contains information such as boundaryconditions, temperature, timesteps, etc. The forcefield is a mathematical functionthat calculates the potential energy of atoms in a system using their location.In the simulation part, force is calculated by differentiating the total energy withrespect to atomic locations, and the atomic locations and velocities are updated.It keeps on calculating forces and updating atomic locations and velocities until amaximum number of timesteps are reached. In the post-processing part, data isanalyzed and then output results are visualized using visualization tools.3Chapter 2Theory and PracticeThis work is based on Datta et al.[8], only the crack geometry is slightly changed andthe set of angles considered. We have used LAMMPS (Large-scale Atomic/MolecularMassively Parallel Simulator) to perform molecular simulations.[9] The dump filesfrom LAMMPS are then used to visualize our results in OVITO.[10]In this work, we have taken a graphene nanoribbon(GNR) of size 150× 150 A° witha pre-crack of half-length 15 A°(a0) and it is tensile loaded in the armchair direction(zigzag graphene).[11] The number of carbon atoms is 8870. It is a 2D simulation,so the z velocity of the atoms is fixed, and to make this a non-periodic simulation,periodic boundary conditions are taken along with a larger size of the simulation boxso that the atoms do not cross the boundary. This is depicted in Figure 2.1.4Figure 2.1: GNR inside the simulation box with a pre-crack oriented at 0°.The x and y axes correspond to the horizontal and vertical axes respectively.2.1 LAMMPSTo perform any MD simulation, typically three files are required. They are inputfiles, data files, and potential files. The data file is visualized in the Figure 2.1 and2.2. Figure 2.2 also shows two orientations of graphene: armchair and zigzag. In thiswork, the GNR is pulled in the armchair direction, so that makes it loaded in zigzag5Figure 2.2: GNR with a pre-crack oriented at 45°.direction, hence it is also called zigzag graphene. We have used AIREBO potentialfor our simulation because it is commonly used for graphene simulations and capturesthe physics better than other potentials available.[12] Metal units are used for all thesimulations. Berendsen thermostat is used to carry out all the simulations at roomtemperature (300K). The timestep given is 0.001 ps and the displacement rate of thetop layer of carbon atoms is 0.01 A°s−1.62.2 Griffith’s brittle fractureZhang, P., Ma, L., Fan, F. et al. have shown that the classic Griffith theory of brittlefracture applies to graphene.[13] The Griffith’s criterion for a central crack of length2a0 is given by,σc =√2γE/πa0 (2.1)where σc is the stress at fracture, E is Young’s modulus, a0 is the half-length of thecrack, and γ is the edge energy. This criterion is used to validate the simulationresults later in this report.7Chapter 3Results and DiscussionThe final image of the MD simulation obtained from OVITO is shown along withthe stress-strain curve in the results below for the corresponding ten angles. Thestress values are then used to calculate the critical stress intensity factor of fracture,Kc. Later, the Kc’s of all ten GNR are then compared with the published results forvalidation.The stress in the images is stress in the Y direction because we are displacing thenanoribbon in the positive y-direction. Rainbow scale is chosen to visualize the stress,where red is the highest, green is in the middle, and blue is the lowest.MATLAB is used to plot all the stress-strain curves.[14] The displacement rate of thetop layer of GNR is taken to be the value that corresponds to clear crack propagation.83.1 GNR with a pre-crack oriented at 0°Figure 3.1: Crack propagation in GNR with a pre-crack at 0°.93.2 GNR with a pre-crack oriented at 5°Figure 3.2: Crack propagation in GNR with a pre-crack at 5°.103.3 GNR with a pre-crack oriented at 10°Figure 3.3: Crack propagation in GNR with a pre-crack at 10°.113.4 GNR with a pre-crack oriented at 15°Figure 3.4: Crack propagation in GNR with a pre-crack at 15°.123.5 GNR with a pre-crack oriented at 20°Figure 3.5: Crack propagation in GNR with a pre-crack at 20°.133.6 GNR with a pre-crack oriented at 25°Figure 3.6: Crack propagation in GNR with a pre-crack at 25°.143.7 GNR with a pre-crack oriented at 30°Figure 3.7: Crack propagation in GNR with a pre-crack at 30°.153.8 GNR with a pre-crack oriented at 35°Figure 3.8: Crack propagation in GNR with a pre-crack at 35°.163.9 GNR with a pre-crack oriented at 40°Figure 3.9: Crack propagation in GNR with a pre-crack at 40°.173.10 GNR with a pre-crack oriented at 45°Figure 3.10: Crack propagation in GNR with a pre-crack at 45°.3.11 Virial stress computation in LAMMPSVirial stress computation as mentioned in Dr. Dibakar Datta’s PhD thesis is givenby equation 3.1.[15] For homogeneous systems, it is a measure of mechanical stress18on molecular level.〈σVij〉=12VN∑k=1N∑l 6=k(xli − xki)fkl( int )j −1VN∑k=1xki(fk(ext)j −mk ¨̂jx)(3.1)The virial stress computation in LAMMPS for an atom numbered 1 is given byequation 3.2.[15] The six components of a symmetric tensor is genereated when a andb take on the values x, y, and z.Sab =−[mvavb +12Np∑n=1(r1aF1b + r2aF2b) +12Nb∑n=1(r1aF1b + r2aF2b) +13Na∑n=1(r1aF1b + r2aF2b + r3aF3b) +14Nd∑n=1(r1aF1b + r2aF2b + r3aF3b + r4aF4b) +14Ni∑n=1(r1aF1b + r2aF2b + r3aF3b + r4aF4b) + Kspace (ria , Fib) +Nf∑n=1riaFib(3.2)This formulation is in pressure(stress) × volume units. To convert it into units ofstress, we need to divide it by volume per-atom. But in a deformed structure, thevolume per-atom is ill defined. So, to get total pressure of the system, the sum ofdiagonal elements of stress tensor of individual atoms for all atoms in the system isdivided by the product of dimension and volume of the system.[15] So, to get stressin GPa, we take a product of LAMMPS value and a multiplication factor.Multiplication factor is given by 1/volume per-atom × 1e-4 (1 bar = 1e-4 GPa). We19have 8860 atoms and size of the nanoribbon is 150 × 150 A°. So, volume per-atomin this work is 8.66 °A3. In our case the multiplication factor is 1.155e-5.[15]3.12 Stress-Strain curveFigure 3.11: Stress-strain curve for all GNRsThe strain values in the Figure 3.11 look higher than expected, but it is dependenton the rate at which the nanoribbon is pulled. In this work, this is checked throughvarying the rate of displacement and observing the change in values of stress at20fracture. No significant change in values of stress at fracture was found, and the rateat which we get a clear crack propagation was taken.3.13 ValidationFrom Zhang, P., Ma, L., Fan, F. et al, the average value of Kc obtained is 4.0 MPa√m,with a standard deviation of 0.6 MPa√m.[13] In Table 3.1, this value is comparedwith Kc values for all ten cracked GNRs. Kc is equal to σc√πa0.[13]Table 3.1Calculated Kc values corresponding to their crack anglescrack angle 0 5 10 15 20 25 30 35 40 45Kc 3.48 3.69 3.45 3.23 3.2 2.9 3.65 3.17 3.59 3.48According to Griffith’s criterion, the product of σc and√a0 has to be constant. σc√a0range in the referred article is 1.73 to 2.78 MPa√m. In this report, only the 25° crackhas a value of 1.64 MPa√m that goes below the given range. Hence, the simulationresults are in accordance with the published results. Table 3.2 contains σc√a0 valuesfor all ten GNRs.Table 3.2Calculated σc√a0 values corresponding to their crack anglescrack angle 0 5 10 15 20 25 30 35 40 45σc√a0 1.96 2.08 1.95 1.82 1.8 1.64 2.06 1.79 2.03 1.9621All the values are around 2 MPa√m, so it can be considered as a constant. Hence,Griffith’s theory is applicable to graphene.3.14 ConclusionAfter analyzing the results from MD simulations, the crack path can be observed forall ten GNRs. It can be concluded that the GNR with crack at 25° angle has thelowest strength because it has the sharpest crack tip and slightly different crack shapethan other GNRs. This conclusion is specific for a particular crack length. It canalso be said that the crack path depends on the length of the crack, orientation ofthe crack, crack tip sharpness, and also the shape of the crack.3.15 Future workSince MD simulations are computationally expensive, our future work involves pre-diction of fracture in graphene using Deep Learning. The goal is to reduce the timerequired in making fracture path predictions by bypassing MD simulations with amachine learning model.22References[1] Wikimedia. AlexanderAlUS. 2010.[2] Büyüköztürk, O.; Buehler, M. J.; Lau, D.; Tuakta, C. International Journal ofSolids and Structures 2011, 48(14), 2131–2140.[3] Lee, C.; Wei, X.; Kysar, J.; Hone, J. Science 2008, 321, 385–388.[4] Balandin, A.A.and Ghosh, S.; Bao, W.; Calizo, I.; Teweldebrhan, D.; Miao, F.;Lau, C. Nano Letters 2008, 8, 902–907.[5] Ho, Y. H.; Wu, J. Y.; Chiu, Y. H.; Wang, J.; Lin, M. F. Royal Society 2010,368, 5353–5354.[6] Wahl, L.; Maas, S.; Waldmann, D.; Zürbes, A.; Frères, P. Journal of SandwichStructures & Materials 2012, 14, 449–468.[7] Molecular dynamics — Wikipedia, the free encyclopedia. Wikipedia contributors.2021.23[8] Datta, D.; Nadimpalli, S. P.; Li, Y.; Shenoy, V. B. Extreme Mechanics Letters2015, 5, 10–17.[9] Plimpton, S. Journal of Computational Physics 1995, 117, 1–19.[10] Stukowski, A. JAN 2010, 18(1).[11] Pei, Q.; Zhang, Y.; Shenoy, V. Carbon 2010, 48(3), 898–904.[12] Stuart, S. J.; Tutein, A. B.; Harrison, J. A. The Journal of Chemical Physics2000, 112(14), 6472–6486.[13] Zhang, P.; Ma, L.; Fan, F.; Zeng, Z.; Peng, C.; Loya, P. E.; Liu, Z.; Gong, Y.;Zhang, J.; Zhang, X.; Ajayan, P. M.; Zhu, T.; Lou, J. Nat Commun 2014, 5,3782.[14] MATLAB. version 9.8.0.1538580 (R2020a); The MathWorks Inc.: Natick, Mas-sachusetts, 2020.[15] Datta, D. Topics in Mechanics at the Bottom, Energy StorageSystems, andEmerging Nanomaterials PhD thesis, Brown University, 2015.24

STUDYING THE EFFECTS OF INITIAL CRACK ANGLE ON THE CRACK PROPAGATION IN GRAPHENE NANO-RIBBON THROUGH MOLECULAR DYNAMICS SIMULATIONS

https://digitalcommons.mtu.edu/cgi/viewcontent.cgi?article=2331&context=etdr

STUDYING THE EFFECTS OF INITIAL CRACK ANGLE ON THE CRACK PROPAGATION IN GRAPHENE NANO-RIBBON THROUGH MOLECULAR DYNAMICS SIMULATIONS

Abstract

Similar works

Full text

Available Versions

Michigan Technological University