Dear lammps users,
Recently, I repeated the lammps script for calculating elastic constant matrix for silicon with sw potential in examples/ELASTIC.
To verify the variable ‘up’ in init.mod, I found that if the origin minimization parameters were used, the calculated results depend on ‘up’.
Define minimization parameters
variable etol equal 0.0
variable ftol equal 1.0e-10
variable maxiter equal 100
variable maxeval equal 1000
variable dmax equal 1.0e-2
For example, I let up = 1e-6,
Elastic Constant C11all = 99.7221839332975 GPa
Elastic Constant C22all = 99.7243000082318 GPa
Elastic Constant C33all = 99.7231207303535 GPa
Elastic Constant C12all = 99.7232385757629 GPa
Elastic Constant C13all = 99.7226489328092 GPa
Elastic Constant C23all = 99.7237069678318 GPa
Elastic Constant C44all = 3.1986914724038e-09 GPa
Elastic Constant C55all = 1.1371837593108e-08 GPa
Elastic Constant C66all = 4.95788640651672e-09 GPa
Elastic Constant C14all = 6.1373976214534e-05 GPa
Elastic Constant C15all = -5.141205612193e-05 GPa
Elastic Constant C16all = 1.57980029290368e-05 GPa
Elastic Constant C24all = 5.32673188937985e-05 GPa
Elastic Constant C25all = -5.372387798094e-05 GPa
Elastic Constant C26all = 1.00255700366896e-05 GPa
Elastic Constant C34all = 5.28559240829355e-05 GPa
Elastic Constant C35all = -5.58359246495825e-05 GPa
Elastic Constant C36all = 1.5894983856242e-05 GPa
Elastic Constant C45all = -1.39644807646324e-08 GPa
Elastic Constant C46all = -1.11850910755457e-09 GPa
Elastic Constant C56all = -9.44928038364181e-09 GPa
if up= 1e-5
Elastic Constant C11all = 93.5568027510024 GPa
Elastic Constant C22all = 93.5098837659731 GPa
Elastic Constant C33all = 93.5314988938877 GPa
Elastic Constant C12all = 93.5331855133706 GPa
Elastic Constant C13all = 93.543993087387 GPa
Elastic Constant C23all = 93.5205338786411 GPa
Elastic Constant C44all = -4.32393211215388e-08 GPa
Elastic Constant C55all = -4.05308377918188e-08 GPa
Elastic Constant C66all = 1.31928797457765e-08 GPa
Elastic Constant C14all = 0.00129211447419147 GPa
Elastic Constant C15all = -0.000409460105511037 GPa
Elastic Constant C16all = 6.63365327868965e-05 GPa
Elastic Constant C24all = 0.00112657241166843 GPa
Elastic Constant C25all = -0.000519625106828757 GPa
Elastic Constant C26all = 0.000267206396890855 GPa
Elastic Constant C34all = 0.00128292271307588 GPa
Elastic Constant C35all = -0.000506465325361475 GPa
Elastic Constant C36all = 3.50020853263245e-05 GPa
Elastic Constant C45all = -1.09219477329496e-09 GPa
Elastic Constant C46all = -1.69407108193573e-10 GPa
Elastic Constant C56all = 8.85300901261172e-10 GPa
if up=1e-4
Elastic Constant C11all = 41.7035734293569 GPa
Elastic Constant C22all = 41.7285159054187 GPa
Elastic Constant C33all = 41.7117441049898 GPa
Elastic Constant C12all = 41.7041160935199 GPa
Elastic Constant C13all = 41.6957160711382 GPa
Elastic Constant C23all = 41.7082092700946 GPa
Elastic Constant C44all = 2.06940578876834e-06 GPa
Elastic Constant C55all = 2.31423174286954e-06 GPa
Elastic Constant C66all = 1.92722081809188e-06 GPa
Elastic Constant C14all = 0.00988657202735753 GPa
Elastic Constant C15all = 0.00345563221637275 GPa
Elastic Constant C16all = -0.00652841666893143 GPa
Elastic Constant C24all = 0.0100700456920097 GPa
Elastic Constant C25all = 0.00363852451487783 GPa
Elastic Constant C26all = -0.00698332973705498 GPa
Elastic Constant C34all = 0.0102229624120423 GPa
Elastic Constant C35all = 0.00372358494640636 GPa
Elastic Constant C36all = -0.00648279202014802 GPa
Elastic Constant C45all = -9.348961524741e-09 GPa
Elastic Constant C46all = 1.48726766543841e-07 GPa
Elastic Constant C56all = 8.93584328939922e-08 GPa
I change the minimization parameters as follows,
Define minimization parameters
variable etol equal 0.0
variable ftol equal 1.0e-10
variable maxiter equal 1000
variable maxeval equal 100000
variable dmax equal 1.0e-2
if let up=1e-4 or 1e-6, the elastic constants are insensitivie to the ‘up’,
Elastic Constant C11all = 50.7116289450506 GPa
Elastic Constant C22all = 50.7116287220582 GPa
Elastic Constant C33all = 50.7116288803589 GPa
Elastic Constant C12all = 50.7014880292042 GPa
Elastic Constant C13all = 50.7014881083634 GPa
Elastic Constant C23all = 50.7014879968484 GPa
Elastic Constant C44all = 2.02880179583673e-06 GPa
Elastic Constant C55all = 2.02906528957374e-06 GPa
Elastic Constant C66all = 2.02867122722663e-06 GPa
Elastic Constant C14all = 1.13527006433756e-10 GPa
Elastic Constant C15all = -1.86670963535907e-10 GPa
Elastic Constant C16all = -1.34190087989925e-11 GPa
Elastic Constant C24all = -8.62558888590386e-11 GPa
Elastic Constant C25all = -1.70246695005236e-10 GPa
Elastic Constant C26all = 4.3461786649661e-11 GPa
Elastic Constant C34all = 8.04265429639297e-11 GPa
Elastic Constant C35all = -7.28952312878578e-11 GPa
Elastic Constant C36all = 4.98277759297564e-11 GPa
Elastic Constant C45all = 2.3365028501589e-11 GPa
Elastic Constant C46all = 1.21848596355367e-10 GPa
Elastic Constant C56all = 3.48368936525295e-11 GPa
But we can find that, all the calculated elastic constants are not corrected.
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
Could anyone point out my mistakes when runing this script. I only modified the init.mod.
variable up equal 1.0e-4
units metal
variable cfac equal 1.0e-4
variable cunits string GPa
variable etol equal 0.0
variable ftol equal 1.0e-10
variable maxiter equal 1000
variable maxeval equal 100000
variable dmax equal 1.0e-2
generate the box and atom positions using a diamond lattice
variable a equal 5.43
boundary p p p
lattice diamond $a
region box prism 0 10.0 0 10.0 0 10.0 0.0 0.0 0.0
create_box 1 box
create_atoms 1 box
Need to set mass to something, just to satisfy LAMMPS
mass 1 28.0855
Thanks very much.
Best Regards,
Mike jiang
June 13, 2018
Prof. Dr Wugui Jiang
M408
School of Aeronautical Manufacturing Engineering
Nanchang Hangkong University