View source for Multi-Scale Modeling of Pure Vanadium

=Abstract=
In this study, the stress/strain behavior of pure vanadium in the elastic and plastic regimes are simulated using a multiscale approach. Density Functional Theory (DFT) calculations were performed to obtain the energy verse volume curves for BCC, FCC, and HCP Vanadium, along with the coherence energy, bulk modulus, generalized stacking fault energy (GSFE), and lattice constants for the stable BCC V phase. Three second nearest-neighbor (2NN) Modified Embedded Atom Method (MEAM) interatomic potential was generated from these results to quantify uncertainty in GSFE result. The MEAM potentials were used in a molecular dynamics code to determine the dislocation mobility. The calculated dislocation mobility was upscaled to Dislocation Dynamics (DD) scale to study the stress-strain response of the material, specifically the hardening behavior. From the DD calculations, hardening law parameters were upscaled to crystal plasticity scale to simulate the stress/strain response of a single and polycrystalline material. Finally, the stress/strain curves calculated from crystal plasticity are used for the parameterization of the internal state variable (ISV) model to enable simulation on a continuum scale.

Authors: Caleb O. Yenusah <sup>a,b</sup>, Abhijith Madabhushi <sup>a</sup>, William M. Furr <sup>a,b</sup>, Benmbarek Mouhsine <sup>a</sup>
<br>
a.	Department of Mechanical Engineering, Mississippi State University, Mississippi State, MS, USA
<br>
b.	Center for Advanced Vehicular Systems, Mississippi State University, Mississippi State, MS, USA

=Introduction=
Integrated Computational Materials Engineering (ICME) sets a framework for the modeling of materials from the lowest length and time scales to the continuum scale. Bridges between the different scales are established in relation to the nature of the problem, or a continuum scale phenomenon of interest. In this case study, the plasticity of Vanadium is the continuum scale phenomenon of interest. To achieve this goal, Density Functional Theory (DFT) calculations were performed in Quantum Espresso<ref>P. Giannozzi et al., “QUANTUM ESPRESSO: a modular and open-source software project for quantum simulations of materials,” J. Phys. Condens. Matter, vol. 21, no. 39, p. 395502, Sep. 2009.</ref> to obtain the energy verse volume curves for BCC, FCC, and HCP Vanadium, along with the coherence energy, elastic constants, generalized stacking fault energy (GSFE), and lattice constants for the stable BCC V phase. To upscale the information generated from the electronic scale to the atomic scale, A second nearest-neighbor (2NN) Modified Embedded Atom Method (MEAM) interatomic potential was generated from these results <ref name='MEAM'>Lee, Byeong-Joo, et al, "Second nearest-neighbor modified embedded atom method potentials for bcc transition metals," Physical Review B, vol. 64, no. 18, 2001.</ref>. This MEAM potential is used in the atomic scale in molecular dynamics (MD) calculations. The MD calculations can then be used to examine atomic scale phenomena like the movement of dislocations as a function of applied stress, temperature, and other boundary conditions. From the MD simulations, the atomic positions in the lattice under the influence of the dislocation are obtained, from which important variables like the dislocation velocities can be calculated; which in turn can be used to obtain the dislocation mobility. The dislocation mobility will be upscaled to the Mesoscale via Dislocation Dynamics (DD) calculations to study the stress-strain response of the material, especially the hardening behavior. From the DD calculations, hardening law parameters were upscaled to crystal plasticity calculations to simulate the stress/strain response of a single and polycrystalline material. Finally, the stress/strain curves calculated from crystal plasticity are used for the parameterization of the internal state variable (ISV) model to enable simulation on a continuum scale. 

This exercise is of utmost importance for capturing the critical stress-strain behavior post yielding known as plastic deformation. During the deformation of a material in the plastic regime, dislocation density increases, resulting in the interaction of dislocation, this interaction dictates the hardening behavior of a material. In this study, the Voce hardening rule parameters for a single crystal are determined from DD calculations, which are then used in a Finite Element Code (FEM) called CPFEM<ref>F. Roters, "Advanced material models for the crystal plasticity finite element method: development of a general CPFEM framework.No. RWTH-CONV-144,865th," in Material Science and Engineering, AACHEN, Germany, 2011</ref> to simulate hardening in a single crystal, followed by polycrystalline simulations for three different loading conditions (Tension, Compression, and Torsion) for a realistic quantification of the stress-strain response of the material to different stress states. This stress-strain behavior of the polycrystalline solid is then used to obtain a parameterized internal state variable model (ISV) for the material. The ISV model serves as a bridge to the continuum scale where large systems of engineering significance can be simulated. 
During the upscaling process, care is taken to capture the sensitivity of the computed values of the current scale to the values from the previous scale. To accomplish this, a sensitivity analysis is performed by creating three different test cases, starting from the MEAM potential creation. The best fit case is the resultant MEAM potential previously mentioned using the MPC calibration tool<ref name='MPC'>Carino, Ricolindo L., and Mark F. Horstemeyer. "Case Studies in Using MATLAB to Build Model Calibration Tools for Multiscale Modeling." Applications from Engineering with MATLAB Concepts. IntechOpen, 2016.</ref>. Information from an upper and lower bound potential is also upscaled to the MD simulations to observe the sensitivity. Similar care is taken at the mesoscale to capture the sensitivity of the values from the atomic scale. This is done to capture the uncertainties in the analyses like the generalized stacking fault energy; which is propagated from the nanoscale to the MD scale calculations of dislocation mobility, which is then upscaled to the mesoscale DD calculations of the hardening behavior, the CPFEM simulation of stress/strain response of single and polycrystalline, and finally, the parameterization of the ISV model. 

While the aim of this research endeavor is to obtain a working knowledge of the ICME design paradigm, future work will involve the improvement of the results produced in this work to accurately capture the continuum scale plasticity behavior of pure BCC Vanadium. The following sections describe the methodology involved in this work and the corresponding results obtained for the best fit, upper bound and lower bound cases.

=Methodology=                                                                                                                                                             <!-- Methodology Section -->
==Electronic Scale Simulations==                                                                                                                                          <!-- Methodology from Homework #1 -->
===Ab Initio calculations===
Ab initio calculations in this work were based on density functional theory (DFT) as implemented in the Quantum Espresso code. Projector augmented-wave (PAW) potential<ref>P. E. Blöchl, “Projector augmented-wave method,” Phys. Rev. B, vol. 50, no. 24, pp. 17953–17979, Dec. 1994.</ref> and the Perdew-Burke-Ernzerhof generalized gradient approximation<ref>J. P. Perdew, K. Burke, and M. Ernzerhof, “Generalized Gradient Approximation Made Simple,” Phys. Rev. Lett., vol. 77, no. 18, pp. 3865–3868, Oct. 1996.</ref> were used. Convergence study was performed in two stages. First to determine the appropriate k-point mesh point needed, the plane wave energy cutoff was fixed at 653.07 eV and k-point mesh were varied. After the converged k-point mesh with respect to the energy cutoff is obtained, a second study was performed using the converged k-point mesh and different energy cutoff ensuring that all results are converged with respect to both parameters. K-point convergence study was also performed for lattice constant and bulk modulus. The results of this study can be found in the Appendix. Energy versus atomic spacing for bcc, fcc, and hcp structures was calculated using the converged energy cutoff and K-point mesh. This is needed for calibration of the MEAM parameters. The lattice constant, bulk modulus, and cohesive energy for the ground state bcc, V was calculated using Murnaghan equation of state (EOS). 
The generalized stacking fault energy (GSFE) values (E<sub>sf</sub>) were calculated using Eq. 1 for calibration of the MEAM parameters for a specific slip system. This ensures that plastic deformation is correctly reproduced by the MEAM potential developed. Due to the fact that vanadium belongs to bcc metal, the (110) plane is the most close-packed plane and the {110}<111> is the most preferable slip system family and should have the smallest GSFE curve. The GSFE curves for (110)[-111] , (110)[001] , and (110)[-110] were plotted. Since the former slip system had the smallest GSFE curve, the MEAM parameters were calibrated to the (110)[-111] slip plane and direction. 
<br>
[[Image:GSFE EQ.JPG|center|250px]]
where E<sub>tot</sub> is the total energy of the structure with a stacking fault, N is the number of atoms in the system, ε is the total energy per atom in the bulk, and A is the unit cell area that is perpendicular to the stacking fault<ref>M. F. Horstemeyer, Integrated Computational Materials Engineering (ICME) for Metals: Using Multiscale Modeling to Invigorate Engineering Design with Science. John Wiley & Sons, 2012.</ref>.

===MEAM parameter calibration===
MEAM parameter calibration was performed using the MEAM Parameter Calibration software (MPC) developed by Mississippi State University<ref name='MPC'></ref>. MEAM parameters were calibrated to the energy verse atomic spacing and GSFE data calculated from DFT and the elastic moduli of bcc V found in the literature<ref>Alers, G. A. "Elastic moduli of vanadium." Physical Review 119, no. 5 (1960): 1532.</ref>. A comparison was made between the DFT calculated energy verse atomic spacing and GSFE curve and results from Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) utilizing the developed MEAM potential.

==Molecular Dynamics Simulations==                                                                                                                                                    <!-- Methodology form Homework #2 -->
The MD simulations were performed using Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS). A “pad” type structure with BCC Vanadium (V) atoms was created for the purpose of the simulations as shown in the figure below. It was periodic in the X and Z, and non-periodic/shrink-wrapped in the Y direction. An Edge dislocation was introduced in the middle of the cell and constant shear stress was applied on the bottom and top face to propagate the dislocation in the X direction. A constant number of atoms, volume, and energy conditions were prescribed. The temperature of the system was fixed at 300K and equilibrated for 10000ps. The dislocation mobility for the edge dislocation as a function of applied stress and velocity was calculated from the MD simulation results. A convergence study was performed for applied shear stress for the different total number of atoms. The results from the MD simulations serve as the bridging parameters to the mesoscale DD calculations to capture the stress-strain behavior of the material. <br>
<br>
[[Image:OVITO.jpg|center|700px]]
<div align="center">Figure 1. Atom positions within the Pad structure after the equilibration stage. </div>

<br>

==MDDP Simulations==                                                                                                                                                      <!-- Methodology form Homework #2 & Homework #3 -->
In the next step, Multiscale Dislocation Dynamics Plasticity (MDDP)<ref>Zbib, Hussein M., Tomas Diaz de la Rubia, "A multiscale model of plasticity," International Journal of Plasticity, vol. 18, no. 9, pp. 1133-1163, 2002</ref> code was used to create a single Frank-Read source (SFRS) to perform the DD calculations. The code used two input files and an executable (BCCdata.exe) to generate the dislocation input file (DDinput). This input file was modified so that the Frank-Read pair created earlier is modified to a single Frank-Read source. Further modifications to the input file to account for the nearest neighbor atom information was corrected. The results of the MDDP calculations were later extracted and visualized using Tecplot. The process was repeated for a single crystal with multiple Frank-Read sources (MFRS). The material properties used in MDDP calculations are detailed in Table 1. Annealed initial state was assumed for the material with an initial dislocation density in the order of 10<sup>12</sup>. Preliminary studies were carried out to understand the effect of the number of Frank-Read sources on the initial dislocation density to determine the number of Frank-Read sources needed for a dislocation density in the order of 10<sup>12</sup>. The MFRS were randomly generated on different planes. For BCC vanadium, all glide planes are of the {110} family. A strain rate of 1.0E5  was used for all simulations. The material properties used in MDDP calculations are detailed in Table 1.
<div align="center">Table 1: Material properties for BCC Vanadium used in MDDP Simulations </div>
[[Image:Matl Props Table.PNG|center|800px]]

All the values in the above table were obtained from the results of the calculations from the previous scales. The Voce Hardening law parameters were obtained from the results of the MDDP simulations which are upscaled to the mesoscale via the CPFEM simulations.

== CPFEM Simulations==                                                                                                                                                    <!-- Methodology from Homework #3 & Homework #4 -->
Crystal plasticity finite element method (CPFEM) simulations were run for single crystals and 500 crystals in tension, compression, and shear. The Voce hardening parameters obtained from the MDDP results were used in these simulations, as well as the parameters in Table 2. C11, C12, and C44 are the single crystal elastic constants obtained from previous MPC calibrations, B is the dislocation drag coefficient previously obtained from Lammps simulations.
<div align="center">Table 2: Material properties for BCC Vanadium used in MDDP Simulations </div>
[[Image:Elastic constants.JPG|center|600px]]
These values were used in an Abaqus user material subroutine (UMAT) and applied to a single element in tension, compression, and shear. The slip mode examined was <-111> {110} for bcc Vanadium. To evaluate the effects of varying the number of crystals, the best fit case obtained from MDDP simulations was simulated in tension and compared.

== Macroscale ISV Model Calibration==                                                                                                                                     <!-- Methodology form Homework #4 -->
The results of the polycrystalline CPFEM simulations were used to generate the ISV model using the MSU ISV DMGfit 55p v1p1 tool<ref>Mississippi State University, "MSU ICME," Mississippi State University, [Online]. Available: https://icme.hpc.msstate.edu/mediawiki/index.php/DMGfit_55p_v1p1. [Accessed Monday April 2019].</ref>. This was done for the purpose of upscaling the results of the mesoscale onto the macroscale simulations. DMGfit was designed to be an interactive calibration tool to model damage and plasticity using the MSU ISV model. The model is a FORTRAN subroutine implemented in ABAQUS as a UMAT (user material subroutine). The purpose of the calibration tool is to determine ‘material constants’ that will be used as inputs in the finite element simulations. The DMG model itself is specified by 55 constants, which include the dataset and the fitted parameters as shown in the figure below. The upscaled results like bulk and shear modulus and other material properties serve as the input to the calibration tool along with the stress-strain data from CPFEM. The tool uses these inputs and subjects them to an optimization protocol called fminsearch. The relevant Best Fit, Upper Bound and the Lower Bound results of the previous scales like the tension, compression and torsion simulations from the polycrystalline CPFEM simulations are saved as data set files and are subjected to this process and three ISV models are obtained. 
The fitted ISV constants are added to an input deck to be used alongside the DMG UMAT to run a single element calculation in ABAQUS for the three loading conditions. The results obtained from the finite element calculation serve as verification of the plasticity-damage calibration.
<br>
[[Image:DMGfit.jpg|center|800px]]
<div align="center">Figure 2. DMGfit 55 v1p1 calibration tool used to fit stress-strain behavior from CPFEM for BCC vanadium </div>

=Results=                                                                                                                                                                 

==Electronic Scale Simulations==
The Shear Modulus, the lattice constant, and the cohesive energy are obtained after the convergence study using the Universal Equations of State (UEOS). Out of the four EOS choices given, BIRCH1 and Keane models were not able to generate results. It seems like the code was unable to converge to a result. From the results generated by BIRCH1 and MURNAGHAN as shown in Table 3, the former yielded a lower lattice constant value but a higher bulk modulus and cohesive energy values. The value generated by MURNAGHAN was closer to Vanadium’s lattice constant as per literature<ref>Molecular dynamics simulation of vanadium using an interatomic potential fitted to finite temperature properties - ScienceDirect.” [Online]. Available: https://www.sciencedirect.com/science/article/pii/S002231150201019X. [Accessed: 09-Mar-2019]</ref>  Therefore MURNAGHAN was used for the rest of the analysis performed in this study.  
<br>  
                                                                                                                                                                                            <!-- Results from Homework #1 -->      
<div align="center">Table 3. EOS study using the converged K-point and energy cut-off. </div>                                                                                                                                   
[[Image:DFT Table1.PNG|center|600px]]

==MEAM Parameter Calibraion==    
As a prelude to this study, the MEAM potential previously obtained was fine-tuned to accurately capture the behavior of the generalized stacking fault energy (GSFE). The coherence energy, bulk modulus, and lattice constants of the stable BCC V phase were well established in the literature from experiments and first principle calculations. However, the characteristics of the GSFE can have only been reported through first principle calculations <ref>Gui, Li-Jiang, and Yue-Lin Liu, "First-principles studying the properties of oxygen in vanadium: Thermodynamics and tensile/shear behavior," Computational Condensed Matter, vol. 7, pp. 7-13, 2016</ref><ref>Zhang, Xingming, et al., "The effects of interstitial impurities on the mechanical properties of vanadium alloys: A first-principles study," Journal of Alloys and Compounds, vol. 701, pp. 975-980, 2017</ref>. Therefore, any errors in it are quantified through a sensitivity analysis, since the GSFE is paramount to the plastic behavior of materials.  This was performed using three different MEAM potentials, a good fit to the DFT calculations, another that underestimated the GSFE curve by 26% <ref name='MEAM'></ref>, and a third one that overestimated the GSFE curve by 22% (Figure 3). These MEAM potentials will hereon be referred to as Best fit, lower bound, and upper bound. It should be noted that the MEAM potential for V from Lee et al [1] underestimated the GSFE curve and this was use used as the lower bound. All the MEAM potential values are shown in Table 4. The values modified in the best fit potential in comparison with the original Lee et al. are A, β(0), β(2), β(3), t(1), t(2), t(3), the attraction and the repulsion parameters. <ref name='MEAM'></ref> MEAM potential to better fit the GSFE curve was therefore calculated from DFT. The energy versus atomic spacing for the three MEAM potentials is reported in Figure A1 in the Appendix along with predicted material properties which were compared with experimental results (Table A1).
<div align="center">Table 4. Parameters for the second nearest-neighbor MEAM potential for V showing upper, lower bounds, and best fit to DFT. </div>                                                                                                                                           
[[Image:MEAM Image1.PNG|center|800px]]
<br>
[[Image:DFT Image1.png|center|500px]]
<div align="center">Figure 3. GSFE Curves for upper, lower bounds, and best fit MEAM to DFT for the [111] direction BCC Vanadium. </div>

==LAMMPS Simulations==                                                                                                                                                      <!-- Results form Homework #2 -->
For the created structure type described in the earlier section, LAMMPS simulations were run to evaluate the dislocation velocity. The LAMMPS output files were then subjected to post-processing using the OVITO visualization software <ref>Alexander Stukowski, "Visualization and analysis of atomistic simulation data with OVITO–the Open Visualization Tool.," Modelling and Simulation in Materials Science and Engineering, vol. 18, no. 1, 2009. </ref>. The OVITO simulations give the position and time values for the edge dislocation during the simulation. The dislocation velocity is obtained by plotting the position value against time and taking its slope.  To select appropriate dimensions for the Pad structure, a convergence study was run to see the behavior of dislocation velocity for applied shear stress of 1500 bar as a function of the number of atoms in the simulation cell. This study gave an optimal number of atoms to run the rest of the simulations. The results of the convergence study and the position vs time plot for the optimal number of atoms (15708) are shown in Figure A3 in the appendix. The value was selected with a view to balancing computational time and the accuracy of the results. A 1-5% range for convergence was considered during the study and the search for convergence was concluded after the last value for the number of atoms shown in Table 3 since a convergence of .457% was attained.
<br>
A graph showing the variation in dislocation velocity with respect to the change in the applied shear stress was plotted (Figure 4) and the characteristics were compared with a similar graph from the ICME for Metals textbook (page no: 388, figure (9.7 a))<ref name='ICME'>Horstemeyer, M. F. (2012). Integrated Computational Materials Engineering (ICME) for metals: using multiscale modeling to invigorate engineering design with science. John Wiley & Sons.</ref>. The figure in the textbook describes the similar effects for an edge location within FCC Aluminum. The velocity scale between the two curves is almost the same however, the velocity of the dislocation in BCC V is lower than FCC Al for the same values of applied stress. Thus, the slope of the FCC Al curve is higher than BCC V. Both curves predict a linear behavior for low applied shear values; however, the linear region in BCC V is between 0-150 MPa which higher than the FCC Al for which the linear region is between 0 – 100 MPa as shown in Figure 4. 
<br> 
[[Image:LAMMPS Image1.png|center|700px]]
<div align="center">Figure 4: (a) Dislocation Velocity vs applied shear stress for edge dislocation in Vanadium vs (b) Figure 9.7a same figure for Al from ICME Textbook </div>
<br>
===Calculation of Dislocation Mobility===
The dislocation mobility of a metal is a very important property that affects its strength and ductility<ref>Kang, Keonwook, Vasily V. Bulatov, and Wei Cai, "Singular orientations and faceted motion of dislocations in body-centered cubic crystals," Proceedings of the National Academy of Sciences, vol. 109, no. 38, pp. 15174-15178, 2012. </ref>. The drag-coefficient (B) is a parameter used to quantify the resistance to dislocation motion. Both quantities were evaluated from the velocity vs applied shear stress graph. The drag-coefficient (B) was found by taking the inverse of the slope of the curve in Figure 4 and multiplying it with the Burgers vector of the edge dislocation. Dislocation mobility, in turn, is the inverse of B. The obtained B and dislocation mobility values are shown in Table 5. 
===Sensitivity Analysis===
After the evaluation of the dislocation mobility and its drag coefficient (B), a sensitivity analysis was conducted on the previously calibrated MEAM potential by comparing it with two other MEAM potentials. The lower bound potential is generated using the MEAM potential for Vanadium evaluated by Lee et al.<ref name='MEAM'></ref> It underestimated the maximum GSFE by 26 %. The upper bound was created using the existing Best fit potential and offsetting its maximum GSFE by 22% more than the best fit. The purpose of the sensitivity analysis is to also study the effect of the GSFE on the Mobility of the Edge dislocation. The LAMMPS simulations were performed for the other two MEAM potentials (upper and lower) using 150 MPa applied shear stress value. Figure 5 shows the position vs time curve for all the three cases and their corresponding velocity vs time curves. The slope for the lower bound potential in the position vs time curve was found to be greater than the best fit and the upper bound curve respectively, indicating higher dislocation velocities for the former.   
[[Image:LAMMPS Image2.png|center|700px]]
<div align="center">Figure 5: (a) Position Vs time curve for Best Fit, Lower Bound and Upper Bound MEAM Potentials atoms, (b) Velocity vs Applied shear stress curve for Best fit, Lower Bound and Upper Bound MEAM potentials (applied shear = 150 MPa)</div>
<br>
<div align="center">Table 5. Dislocation Drag Coefficient and Dislocation Mobilities for the available MEAM Potentials</div>
[[Image:LAMMPS Image3.PNG|center|700px]]
<br>

==MDDP Simulations==                                                                                                                                                        <!-- Results form Homework #2 & Homework #3 -->
Figure 6(a-c) shows the dislocation motion under constant applied stress for the best fit MEAM.  Figure 6a shows the initial SFRS, the propagation of the dislocation as the strain increases is shown in these figures.
[[Image:MDDP Image1.png|center|700px]]
<div align="center">Figure 6: MDDP Simulations for Single FRS for best fit MEAM at different frames (a) 0 (b) 27 (c) 91</div>
<br>
Multiple Frank-Read sources were modeled next to get the dislocation evolution and hardening behavior. Dislocation density evolution for different number of FRS is shown in Figure 7a. The initial dislocation density was 6.30E12, 1.30E13, 1.88E13, and 2.55E13 for 10, 20, 30, and 40 FRS respectively. The 10 FRS was shown for the subsequent calculation for the assumed anneal initial state of the material having dislocation density of the order of 1E12. From Figure 1, forest hardening region was assumed at 7E-8s for 10FRS, while the others have not yet reached a saturation point. The stress/strain response of the single crystal shown in Figure 7b. The stress/strain curve is identical for all the simulations with a yield point between 30MPa to 35MPa. The yield strain was about 3E-2m/m for all initial number of FRS.
[[Image:MDDP Image2.png|center|800px]]
<div align="center">Figure 7: (a) Dislocation density (DD) vs. time for different numbers of Frank-Read sources (b) Stress vs. strain curve for different numbers of Frank-Read sources </div>
===Voce Hardening===
The Voce hardening equation can be written as:
[[Image:Voce1.JPG|center|350px]]
where 𝜅<sub>𝑠</sub>,𝜅<sub>0</sub>,𝑎𝑛𝑑 ℎ<sub>0</sub> are material properties obtained by correlating the equation with the hardening evolution predicted by DD. Because DD calculations start initially using a random distribution of source segments, the dislocation density calculated at the beginning of the simulation can be ignored and only the forest hardening region, which captures dislocation interactions would be used in the determination of the Voce hardening relationship. The dislocation hardening can be written using the classical Taylor relation:
[[Image:Voce2.JPG|center|250px]]
Here, ρ<sub>f</sub> is the forest dislocation density, b is the magnitude of the Burgers vector, μ shear modulus, and α is a constant representing an average of the junction strength over all existing configurations. α=0.35 in this work<ref name='ICME'></ref>. Using Eq. 3, hardening (κ) as a function of time for various initial FRS is calculated and the results are shown in table 6 along with the fitted Voce equation and material parameters κ<sub>s</sub>,κ<sub>0</sub>, and h<sub>0</sub>. A forest hardening region beginning at about 7E-8s can be observed for 10 FRS. For the others, a forest hardening region could not be identified, and extended simulation times are needed to identify the hardening region. But from a rough fit, κ<sub>0</sub> increases with the number of FRS.
<div align="center">Table 6. Fitted Voce Hardening parameters for Best fit, Upper Bound and Lower Bound for 10 FRS.</div>
[[Image:VOCE parameters.JPG|center|600px]]
<br>

== CPFEM Simulations==                                                                                                                                                      <!-- Results from Homework #3 & Homework #4 -->
The single crystal CPFEM simulations ran for all three cases are presented below in Figure A5 in the appendix. All three cases look identical to the same loading condition. This is expected for single crystal simulations, as the effects of dislocation hardening will not be seen for single crystal simulations. Polycrystalline simulations run with 500 grains were also simulated in tension, compression, and shear for all three cases. The stress-strain plots for these simulations are presented in Figure 7. It is apparent that the different dislocation hardening parameters for each case have affected these results since the yield stress obtained from this exercise is significantly different in comparison with the yield stress observed in BCC Vanadium.
[[Image:CPFEM Image1.png|center|800px]]
<div align="center">Figure 7. CPFEM results for 500 crystal simulations in tension, compression, and torsion for (a) the best-fit case, (b) the lower bound, and (c) the upper bound case.</div>

== Macroscale ISV Model Calibration==                                                                                                                                       <!-- Results form Homework #4 -->
The ISV model for the three loading conditions was calibrated using the DMGfit 55 v1 p1 tool. The three datasets for each case are loaded on the calibration tool with the prescribed material constants and the stress-strain information from CPFEM. Since the data is temperature and strain rate independent, the fitting of the parameters is done in the following order: Compression, Tension, and Torsion. In the compression dataset,  the first constant to be fitted is C3 which relates to Yield, followed by the Kinematic hardening and recovery constants C9 and C7, and finally the Isotropic hardening and recovery constants C15 and C13. After this, the Tension dataset is restored and the process is repeated for the Torsion data set. For the torsion data set, the torsion, compression and tension differentiation constants Ca and Cb are adjusted. This is done for the Best fit, Upper bound and the Lower bound cases. The fitted constants are then added to an input deck to be used to run single element calculations in ABAQUS. The following table lists the ISV fitted constants for the best-fit, upper bound and lower bound cases.
<div align="center">Table 6. Verification results where best fit, lower bound and upper bound columns are organized left to right.</div>
[[Image:DMGfittable.JPG|center|800px]]
<br>
Verification simulations ran with the calibrated model in single element simulations are presented in Figure 8 alongside the calibration plots from the DMGfit tool. These plots show reasonable agreement between the two codes, but differences are apparent. This could be due to several different reasons. Apart from the uncertainty in the GSFE from the previous scales, the current ISV model does not consider temperature and strain rate dependencies. Moreover, the accuracy of the Voce Hardening law parameters obtained in the previous section is the primary concern, since it would take a significantly higher amount of computational power and time to capture the entire scope of the MDDP simulations. These issues could be handled for a more robust study which is outside the scope of this report. 
[[Image:ISV Calib Image1.png|center|1000px]]
<div align="center">Figure 8: Verification results where best fit, lower bound, and upper bound columns are organized left to right.</div>

=Conclusion=<!-- Conclusion Section -->
In this study, we have demonstrated a multiscale approach for the study of stress-strain response of a material, using pure vanadium for illustration. Starting from the electronic scale, DFT calculations provided information about energy verse volume curves, coherence energy, bulk modulus, GSFE, and lattice constant. This information was used in the parameterization of three Nearest-neighbor (2NN) MEAM interatomic potential to capture uncertainty in GSFE result. The MEAM potentials were used in MD code (LAMMPS) to determine the dislocation mobility. The calculated dislocation mobility was upscaled to the DD scale to study the stress/strain response of the material, specifically the hardening behavior. From the DD calculations, hardening law parameters were upscaled to crystal plasticity scale to simulate the stress-strain response of a single and polycrystalline material. Finally, the stress-strain curves calculated from crystal plasticity are used for the parameterization of the internal state variable (ISV) model to enable simulation on a continuum scale.
A more robust study can be performed later for more accuracy leading to the improvement of the results produced in this work to accurately capture the continuum scale plasticity behavior of pure BCC Vanadium. The page will be updated in due course of time as we obtain newer and better results that accurately predict the continuum scale stress-strain behavior of pure BCC Vanadium.

=References=<!-- References Section -->
<references/>

=Appendix=                                                                                                                                                                  <!-- Appendix Section -->
[[Image:EVcurve.jpg|center|1000px]]
<br>
<div align="center">Figure A1. DFT results: MEAM Calibration: Energy vs. Atomic spacing for BCC, FCC, and HCP predicted by the three MEAM potentials (a) Upper bound (b) Best Fit (c) Lower bound.
<br>
<br>
<div align="center">Table A1. Comparison of the MEAM potentials to properties of V reported from experiments. All Experimental data are taken from Lee et al.[2]  </div>
[[Image:TABLE A1.JPG|center|700px]]
<br>
[[Image:Appendix1.png|center|700px]]
<br>
<div align="center">Figure A2. DFT results: Energy vs Atomic spacing curves for (a) various K-point meshes for constant energy cut-off value (b) various energy cut-off for constant converged K-point.</div>
<br>
<br>

[[Image:Appendix3.PNG|center]]
<br>
<div align="center">Figure A3. LAMMPS simulation results: Convergence graph showing convergence at N=15708</div>
<br>
<br>

[[Image:Appendix4.png|center|400px]]
<br>
<div align="center">Figure A4. LAMMPS simulation results: Position vs time curves for various applies shear stresses for N =15708</div>
<br>
<br>

[[Image:Appendix5.png|center|700px]]
<br>
<div align="center">Figure A5. CPFEM results for single crystal simulations in tension, compression, and torsion for (a) the best fit case, (b) the lower bound, and (c) the upper bound case.</div>
<br>
<br>

[[Image:Appendix6.PNG|center|500px]]
<br>
<div align="center">Figure A6. CPFEM results for varying numbers of grains in tension with best fit parameters.</div>

[[Category: Research Paper]]
[[Category: macroscale]]
[[Category: mesoscale]]
[[Category: microscale]]
[[Category: nanoscale]]
[[Category: Electronic Scale]]
[[Category: metals]]