in.elastic
b019664c867elammm-devel

in.elastic
View Options

			# Compute elastic constant tensor for a crystal
			#
			# This script uses the following three include files.
			#
			# init.mod (must be modified for different crystal structures)
			# Define units, deformation parameters and initial
			# configuration of the atoms and simulation cell.
			#
			#
			# potential.mod (must be modified for different pair styles)
			# Define pair style and other attributes
			# not stored in restart file
			#
			#
			# displace.mod (displace.mod should not need to be modified)
			# Perform positive and negative box displacements
			# in direction ${dir} and size ${up}.
			# It uses the resultant changes
			# in stress to compute one
			# row of the elastic stiffness tensor
			#
			# Inputs variables:
			# dir = the Voigt deformation component
			# (1,2,3,4,5,6)
			# Global constants:
			# up = the deformation magnitude (strain units)
			# cfac = conversion from LAMMPS pressure units to
			# output units for elastic constants
			#
			#
			# To run this on a different system, it should only be necessary to
			# modify the files init.mod and potential.mod. In order to calculate
			# the elastic constants correctly, care must be taken to specify
			# the correct units in init.mod (units, cfac and cunits). It is also
			# important to verify that the minimization of energy w.r.t atom
			# positions in the deformed cell is fully converged.
			# One indication of this is that the elastic constants are insensitive
			# to the choice of the variable ${up} in init.mod. Another is to check
			# the final max and two-norm forces reported in the log file. If you know
			# that minimization is not required, you can set maxiter = 0.0 in
			# init.mod.
			#

			include init.mod
			include potential.mod

			# Compute initial state
			fix 3 all box/relax aniso 0.0
			minimize ${etol} ${ftol} ${maxiter} ${maxeval}

			variable tmp equal pxx
			variable pxx0 equal ${tmp}
			variable tmp equal pyy
			variable pyy0 equal ${tmp}
			variable tmp equal pzz
			variable pzz0 equal ${tmp}
			variable tmp equal pyz
			variable pyz0 equal ${tmp}
			variable tmp equal pxz
			variable pxz0 equal ${tmp}
			variable tmp equal pxy
			variable pxy0 equal ${tmp}

			variable tmp equal lx
			variable lx0 equal ${tmp}
			variable tmp equal ly
			variable ly0 equal ${tmp}
			variable tmp equal lz
			variable lz0 equal ${tmp}

			# These formulas define the derivatives w.r.t. strain components
			# Constants uses $, variables use v_
			variable d1 equal -(v_pxx1-${pxx0})/(v_delta/v_len0)*${cfac}
			variable d2 equal -(v_pyy1-${pyy0})/(v_delta/v_len0)*${cfac}
			variable d3 equal -(v_pzz1-${pzz0})/(v_delta/v_len0)*${cfac}
			variable d4 equal -(v_pyz1-${pyz0})/(v_delta/v_len0)*${cfac}
			variable d5 equal -(v_pxz1-${pxz0})/(v_delta/v_len0)*${cfac}
			variable d6 equal -(v_pxy1-${pxy0})/(v_delta/v_len0)*${cfac}

			# Write restart
			write_restart restart.equil

			# uxx Perturbation

			variable dir equal 1
			include displace.mod

			# uyy Perturbation

			variable dir equal 2
			include displace.mod

			# uzz Perturbation

			variable dir equal 3
			include displace.mod

			# uyz Perturbation

			variable dir equal 4
			include displace.mod

			# uxz Perturbation

			variable dir equal 5
			include displace.mod

			# uxy Perturbation

			variable dir equal 6
			include displace.mod

			# Output final values

			variable C11all equal ${C11}
			variable C22all equal ${C22}
			variable C33all equal ${C33}

			variable C12all equal 0.5*(${C12}+${C21})
			variable C13all equal 0.5*(${C13}+${C31})
			variable C23all equal 0.5*(${C23}+${C32})

			variable C44all equal ${C44}
			variable C55all equal ${C55}
			variable C66all equal ${C66}

			variable C14all equal 0.5*(${C14}+${C41})
			variable C15all equal 0.5*(${C15}+${C51})
			variable C16all equal 0.5*(${C16}+${C61})

			variable C24all equal 0.5*(${C24}+${C42})
			variable C25all equal 0.5*(${C25}+${C52})
			variable C26all equal 0.5*(${C26}+${C62})

			variable C34all equal 0.5*(${C34}+${C43})
			variable C35all equal 0.5*(${C35}+${C53})
			variable C36all equal 0.5*(${C36}+${C63})

			variable C45all equal 0.5*(${C45}+${C54})
			variable C46all equal 0.5*(${C46}+${C64})
			variable C56all equal 0.5*(${C56}+${C65})

			# For Stillinger-Weber silicon, the analytical results
			# are known to be (E. R. Cowley, 1988):
			# C11 = 151.4 GPa
			# C12 = 76.4 GPa
			# C44 = 56.4 GPa

			print "Elastic Constant C11all = ${C11all} ${cunits}"
			print "Elastic Constant C22all = ${C22all} ${cunits}"
			print "Elastic Constant C33all = ${C33all} ${cunits}"

			print "Elastic Constant C12all = ${C12all} ${cunits}"
			print "Elastic Constant C13all = ${C13all} ${cunits}"
			print "Elastic Constant C23all = ${C23all} ${cunits}"

			print "Elastic Constant C44all = ${C44all} ${cunits}"
			print "Elastic Constant C55all = ${C55all} ${cunits}"
			print "Elastic Constant C66all = ${C66all} ${cunits}"

			print "Elastic Constant C14all = ${C14all} ${cunits}"
			print "Elastic Constant C15all = ${C15all} ${cunits}"
			print "Elastic Constant C16all = ${C16all} ${cunits}"

			print "Elastic Constant C24all = ${C24all} ${cunits}"
			print "Elastic Constant C25all = ${C25all} ${cunits}"
			print "Elastic Constant C26all = ${C26all} ${cunits}"

			print "Elastic Constant C34all = ${C34all} ${cunits}"
			print "Elastic Constant C35all = ${C35all} ${cunits}"
			print "Elastic Constant C36all = ${C36all} ${cunits}"

			print "Elastic Constant C45all = ${C45all} ${cunits}"
			print "Elastic Constant C46all = ${C46all} ${cunits}"
			print "Elastic Constant C56all = ${C56all} ${cunits}"

lammps/examples/elastic/in.elasticb019664c867elammm-devel

in.elasticView Options

lammps/examples/elastic/in.elastic
b019664c867elammm-devel

in.elastic
View Options