CAVS / Capabilities / Materials Modeling / ICME / Code Repository / Magnesium Grain

Integrated Computational Materials Engineering (ICME)

Grain boundary generation in Mg

Abstract

The objective of this research is to generate grain boundary structures over a wide range of degrees of freedom for future use in assessing how grain boundary degrees of freedom impact the properties of polycrystalline Mg. For each grain boundary, this work uses a parallel molecular dynamics code, LAMMPS, with in-plane translations and atom deletion criteria to sample a large number of potential structures to find the global minimum energy grain boundary structure[1][2]. The significance of this research is that grain boundary properties play an important role in the properties of polycrystalline materials and this research enables future atomistic research investigating these properties.

Author(s): Chesley G. Rhodes, Mark A. Tschopp

Corresponding Author: Mark Tschopp

Methodology

The following input script shows how multiple translations and an atom deletion criteria are used to calculate the minimum energy structure. This input script for LAMMPS can be called with a command of the form, "lmp_exe < input.script." This script contains loops over x-translations, z-translations, and atom overlap distances (an atom is deleted when an atom pair with a nearest neighbor distance is less than this distance). If script is saved as a text document, it should be named grain_boundary_simplified.in.

###########################################################################
# __             ___       __        __   __        __   ______    ______ #
#|  |           / _ \     |  \      /  | |  \      /  | |  __  |  |  ____|#
#|  |          / /_\ \    |   \    /   | |   \    /   | | |__|  | | |____ #
#|  |         /  ___  \   |    \  /    | |    \  /    | |   ___|  |____  |# 
#|  |_____   /  /   \  \  |  |\ \/ /|  | |  |\ \/ /|  | |  |       ____| |#
#|________| /__/     \__\ |__| \__/ |__| |__| \__/ |__| |__|      |______|#
###########################################################################

###########################################################################
################SYMMETRIC TILT GRAIN BOUNDARY ENERGY SCRIPT################
###############This script requires an input tilt direction,###############
##############normal direction, and other varibales in order ##############
#############to create a simulation with two grains stacked on#############
############top of each other with a boundary in in the x-z    ############
###########direction. The lower boundary is then moved relative ###########
##########to the upper boundary in the x and z directions, and   ##########
#########delete overlap distance is increased, all to find the    #########
########lowest energy and thus most stable structure. LAMMPS dumps ########
#######the GB Energy for every configuration it tries, with the     #######
######counter number and GBE as the two numbers in the file name to  ######
#####simplify post processing.                                        #####
###########################################################################

###########################################################################
######################Required Input-Orientations##########################
###########################################################################

#Tilt, normal, and periodic directions are required,
#and must be (at least very close to) mutally orthogonal
#in all three dimensions in order to produce good results


#NOTE: This script will only work with the [0 0 0 1] tilt direction, and
#assumes an initial periodic direction of [2 -1 -1 0] (that is, the periodic
#direction if the misorientation angle is equal to 0).

#Tilt direction, in the form of [t1 t2 t3 t4]

variable t1 equal 0
variable t2 equal 0
variable t3 equal -(${t1}+${t2})
variable t4 equal 1

#Normal direction, in the form of [n1 n2 n3 n4]

variable n1 equal 7
variable n2 equal 10
variable n3 equal -(${n1}+${n2})
variable n4 equal 0

#Periodic direction, in the form of [p1 p2 p3 p4]

variable p1 equal 9
variable p2 equal -8
variable p3 equal -(${p1}+${p2})
variable p4 equal 0

###########################################################################
#########################Required Input-Variables##########################
###########################################################################

variable latticeconst equal 3.181269601
variable caratio equal 1.632993162
variable minimumenergy equal -1.528662571
variable xtranslations equal 5
variable ztranslations equal 5
variable pairstyle index "eam/fs"
variable paircoeff index "* * Mg_mm.eam.fs Mg Mg"

#Desired filename for LAMMPS to create and dump atoms to
variable filename index "minimum_energy"

###########################################################################
#########################Direction transformations#########################
###########################################################################

#Convert directions into xyz coordinates
variable f1 equal ${p1}-0.5*${p2}-0.5*${p3}
variable f2 equal sqrt(3)/2*${p2}-sqrt(3)/2*${p3}
variable f3 equal ${p4}*${caratio}

variable g1 equal ${n1}-0.5*${n2}-0.5*${n3}
variable g2 equal sqrt(3)/2*${n2}-sqrt(3)/2*${n3}
variable g3 equal ${n4}*${caratio}

variable h1 equal ${t1}-0.5*${t2}-0.5*${t3}
variable h2 equal sqrt(3)/2*${t2}-sqrt(3)/2*${t3}
variable h3 equal ${t4}*${caratio}

variable f4 equal ${f1}
variable f5 equal -1*${f2}
variable f6 equal ${f3}

variable g4 equal -1*${g1}
variable g5 equal ${g2}
variable g6 equal ${g3}

variable h4 equal ${h1}
variable h5 equal ${h2}
variable h6 equal ${h3}

###########################################################################
####################Calculate a1 a2 a3 components##########################
###########################################################################

#This will be used for inducing the rotation whenever the custom lattice command
#is implemented. It has to be done in this method instead of "orients" due to 
#magnesium's non-repeating, non-terminating c/a ratio causing LAMMPS to crash
#on its othogonality checks. This method is more forgiving, and any round-off
#errors are negated when the initial minimization is run later.

# Calculate lattice parameters  

variable xlattice equal ${latticeconst}
variable ylattice equal sqrt(3)
variable zlattice equal ${caratio}

variable fmag equal sqrt((${f1})^2+(${f2})^2+(${f3})^2)
variable gmag equal sqrt((${g1})^2+(${g2})^2+(${g3})^2)
variable hmag equal sqrt((${h1})^2+(${h2})^2+(${h3})^2)

#For upper Grain Boundary
variable fx equal ${f1}/${fmag}
variable fy equal ${g1}/${gmag}
variable fz equal ${h1}/${hmag}
variable gx equal ${f2}*${ylattice}/${fmag}
variable gy equal ${g2}*${ylattice}/${gmag}
variable gz equal ${h2}*${ylattice}/${hmag}
variable hx equal ${f3}*${zlattice}/${fmag}
variable hy equal ${g3}*${zlattice}/${gmag}
variable hz equal ${h3}*${zlattice}/${hmag}

#For Lower Grain Boundary
variable fx2 equal ${f4}/${fmag}
variable fy2 equal ${g4}/${gmag}
variable fz2 equal ${h1}/${hmag}
variable gx2 equal ${f5}*${ylattice}/${fmag}
variable gy2 equal ${g5}*${ylattice}/${gmag}
variable gz2 equal ${h2}*${ylattice}/${hmag}
variable hx2 equal ${f6}*${zlattice}/${fmag}
variable hy2 equal ${g6}*${zlattice}/${gmag}
variable hz2 equal ${h3}*${zlattice}/${hmag}

###########################################################################
############################Simulation Parameters##########################
###########################################################################

variable etol equal 1.0e-25 
variable ftol equal 1.0e-25 
variable maxiter equal 10000
variable maxeval equal 100000
variable overlapboth equal 2 
variable counter equal 0 

#Calculate simulation box size

variable xlen equal ${xlattice}*sqrt((${f1})^2+(${f2})^2+(${f3})^2)
variable ylen equal ${xlattice}*sqrt((${g1})^2+(${g2})^2+(${g3})^2)
variable zlen equal ${xlattice}*sqrt((${h1})^2+(${h2})^2+(${h3})^2)

#Find y period length to mitigate the affect of the the 
#complimentary GB caused by y periodic conditions
#Simulations run to determine required length converged 
#on 160 nm or greater

variable mult equal "ceil(160/v_ylen)"
variable ylengb equal "v_ylen * v_mult"
variable negylengb equal "-1*v_ylengb"
print "xlen = ${xlen}, ylen = ${ylengb}, zlen = ${zlen}"

# Implement overlap criterion
variable overlapinc equal 531

#Create directory to dump atoms
shell mkdir ${filename}

###########################################################################
####################Begin Loop Structure###################################
###########################################################################
#This will allow multiple configurations as defined by the
#Xtranslations,Ztranslations, and Overlap distance

label loopa 
variable a loop ${xtranslations} 
variable tx equal "(v_a-1)/v_xtranslations*v_xlen" 
    label loopb 
    variable b loop ${ztranslations} 
    variable tz equal "(v_b-1)/v_ztranslations*v_zlen" 
        label loopd 
        variable d loop ${overlapboth} 
            label loopc 
            variable c loop ${overlapinc} 
            variable overlapdist equal ".500 + 0.005 * (v_c-1)" 
 
#Calculate counter for each loop
            variable ctemp equal ${counter}+1 
            variable counter equal ${ctemp} 
            variable ctemp delete 
            print "Counter: ${counter}" 
            
###########################################################################
#####################Initiallize Simulation################################
###########################################################################

            clear
            units metal
            dimension	3
            boundary	p	p	p
            atom_style	atomic
            atom_modify map array

###########################################################################
#######################Create Lattice Structure############################
###########################################################################

            region	box block 0 ${xlen} ${negylengb} ${ylengb} 0 ${zlen} units box
            create_box	2 box
            lattice custom ${xlattice} a1 ${fx} ${fy} ${fz} a2 ${gx} ${gy} ${gz} a3 ${hx} ${hy} ${hz} basis 0.0 0.0 0.0 basis 0.5 0.5 0 basis 0.5 0.83333333 0.5 basis 0 0.33333333 0.5 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
            region	top block INF INF 0 INF INF INF units box
            create_atoms 1 region top basis 1 1 basis 2 1 basis 3 1 basis 4 1  

            lattice custom ${xlattice} a1 ${fx2} ${fy2} ${fz2} a2 ${gx2} ${gy2} ${gz2} a3 ${hx2} ${hy2} ${hz2} basis 0.0 0.0 0.0 basis 0.5 0.5 0 basis 0.5 0.83333333 0.5 basis 0 0.33333333 0.5 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1 
            region  bottom block INF INF INF 0 INF INF units box
            create_atoms 2 region bottom basis 1 2 basis 2 2 basis 3 2 basis 4 2  
            group upper type 1
            group lower type 2

###########################################################################
#########################Induce Pair Potentials############################
###########################################################################

            pair_style	${pairstyle}
            pair_coeff	${paircoeff}
            neighbor 2.0 bin
            neigh_modify delay 10 check yes

###########################################################################
###############Translate atoms and delete overlapping atoms################
###########################################################################

            displace_atoms upper move ${tx} 0 ${tz} units lattice 
            if $c == 1 then "variable atomprev equal 1" 
            if $d == 1 then "delete_atoms overlap ${overlapdist} lower upper"
            if $d == 2 then "delete_atoms overlap ${overlapdist} upper lower" 

#Ensure that equivalent configuration has not been tested
#If the same number of atoms is present as before, keep looping
#until overlap distance increases enough to delete more atoms

            variable natoms equal "count(all)" 
            print "Previous: ${atomprev}, Present: ${natoms}" 
            if ${atomprev} == ${natoms} then "jump grain_boundary_simplified.in loopend" 
 
###########################################################################
#############################Computes######################################
###########################################################################

            compute csym all centro/atom 
            compute eng all pe/atom 
            compute eatoms all reduce sum c_eng 

###########################################################################
########################Initial Minimization###############################
###########################################################################

            reset_timestep 0 
            thermo 100
            thermo_style custom step pe lx ly lz press pxx pyy pzz c_eatoms
            min_style cg
            minimize 1e-15 1e-15 ${maxiter} ${maxeval}
            
 
#Calculate GB Energy after first minimization

            variable gbarea equal "lx * lz * 2" 
            variable gbe equal "(c_eatoms - (v_minimumenergy * count(all)))/v_gbarea" 
            variable gbemJm2 equal ${gbe}*16021.7733 
            variable gbernd equal round(${gbemJm2}) 
            print "After first minimization:" 
            print "GB energy is ${gbemJm2} mJ/m^2" 

###########################################################################
##############################Second Minimization##########################
###########################################################################

#This time allow box bound to expand and relax while minimizing
#and force convergence tolerance to decrease

            reset_timestep 0
            thermo 100 
            thermo_style custom step pe lx ly lz press pxx pyy pzz c_eatoms 
            min_style cg
            fix 1 all box/relax aniso NULL 0.0 NULL vmax 0.0001
            minimize ${etol} ${ftol} ${maxiter} ${maxeval}
 
 
#Calculate GB Energy after first minimization 

            variable gbarea equal "lx * lz * 2" 
            variable gbe equal "(c_eatoms - (v_minimumenergy * count(all)))/v_gbarea" 
            variable gbemJm2 equal ${gbe}*16021.7733 
            variable gbernd equal round(${gbemJm2}) 
            print "After second minimization:" 
            print "GB energy is ${gbemJm2} mJ/m^2" 

# Store number of atoms for overlap criterion for use in comparing next configurations 
            variable atomprev equal ${natoms}

###########################################################################
#####################Dump data into Data file##############################
###########################################################################

            shell cd ${filename}
            reset_timestep 0 
            timestep 0.001 
            dump 1 all custom 5000 dump.${counter}_${gbernd} id type x y z c_csym c_eng 
            run 0             
            shell cd ..

###########################################################################
#####################Finish Loop Structure#################################
###########################################################################

            label loopend 
            next c 
            jump grain_boundary_simplified.in loopc 
        variable c delete 
        next d 
        jump grain_boundary_simplified.in loopd 
    variable d delete 
    next b 
    jump grain_boundary_simplified.in loopb 
variable b delete 
next a 
jump grain_boundary_simplified.in loopa 
print "All done"

Acknowledgments

The authors would like to acknowledge the support from Department of Energy, Southern Regional Center for Innovative Design (SRCLID) program, Contract No.: DE-FC26-06NT42755.

References

The initial methodology was used in the following papers:

  1. ↑ Tschopp, M. A., & McDowell, D.L. (2007). Structures and energies of Sigma3 asymmetric tilt grain boundaries in Cu and Al. Philosophical Magazine, 87, 3147-3173 (http://dx.doi.org/10.1080/14786430701455321).
  2. ↑ Tschopp, M. A., & McDowell, D.L. (2007). Asymmetric tilt grain boundary structure and energy in copper and aluminum. Philosophical Magazine, 87, 3871-3892 (http://dx.doi.org/10.1016/j.commatsci.2010.02.003).