420 likes | 726 Views
South China University of Technology. Structural Relaxation. Xiaobao Yang Department of Physics. http://www.compphys.cn/~xbyang/lectures.html. www.compphys.cn. Schrödinger Equation. INCAR POTCAR KPOINTS POSCAR. KPOINTS Automatic mesh
E N D
South China University of Technology Structural Relaxation Xiaobao Yang Department of Physics http://www.compphys.cn/~xbyang/lectures.html www.compphys.cn
INCAR POTCAR KPOINTS POSCAR KPOINTS Automatic mesh 0 ! number of k-points = 0 ->automatic generation scheme Monkhorst-Pack ! select Monkhorst-Pack (first letter is significant) 1 1 1 ! size of mesh (4x4x4 points along b1, b2, b3) 0. 0. 0. ! shift of the k-mesh INCAR ENCUT = 700 ISMEAR = 0 SIGMA =0.2 NSW = 0 EDIFF = 0.001 POSCAR H2 1.0 12 0 0 0 12 0 0 0 12 2 Cart 0 0 0 0 0 0.75 POTCAR PAW_PBE H 15Jun2001 1.00000000000000000 parameters from PSCTR are: VRHFIN =H: ultrasoft test LEXCH = PE EATOM = 12.4884 eV, .9179 Ry
running on 2 nodes distr: one band on 1 nodes, 2 groups vasp.4.6.31 08Feb07 complex POSCAR found : 1 types and 2 ions LDA part: xc-table for Pade appr. of Perdew POSCAR, INCAR and KPOINTS ok, starting setup WARNING: wrap around errors must be expected FFT: planning ... 1 reading WAVECAR WARNING: random wavefunctions but no delay for mixing, default for NELMDL entering main loop N E dE d eps ncg rms rms(c) DAV: 1 -0.469211785182E+01 -0.46921E+01 -0.23142E+02 12 0.981E+01 DAV: 2 -0.738346878115E+01 -0.26914E+01 -0.26914E+01 12 0.252E+01 DAV: 3 -0.741841946701E+01 -0.34951E-01 -0.34951E-01 18 0.423E+00 DAV: 4 -0.741843794358E+01 -0.18477E-04 -0.18477E-04 12 0.955E-02 DAV: 5 -0.741843809963E+01 -0.15606E-06 -0.15606E-06 18 0.878E-03 0.436E+00 DAV: 6 -0.698113710837E+01 0.43730E+00 -0.22006E-01 12 0.179E+00 0.237E+00 DAV: 7 -0.674222335613E+01 0.23891E+00 -0.52473E-01 24 0.281E+00 0.269E-01 DAV: 8 -0.674568380323E+01 -0.34604E-02 -0.54198E-03 12 0.257E-01 0.257E-01 DAV: 9 -0.674122836271E+01 0.44554E-02 -0.74399E-03 12 0.325E-01 0.324E-02 DAV: 10 -0.674397109933E+01 -0.27427E-02 -0.55758E-04 12 0.997E-02 0.187E-02 DAV: 11 -0.674418356913E+01 -0.21247E-03 -0.58782E-05 18 0.303E-02 1 F= -.67441836E+01 E0= -.67441836E+01 d E =-.241108E-10 writing wavefunctions
Test everything to ensure convergence cluster, nanowire, thin film
for j in 250 300 350 400 450 500 600 700 do cat >INCAR<<! ENCUT = $j ISMEAR = 0 SIGMA =0.2 NSW = 0 EDIFF = 0.001 ! for i in 60 65 70 75 80 85 90 do cat >POSCAR<<! H2 1.0 12 0 0 0 12 0 0 0 12 2 Cart 0 0 0 0 0 0.$i ! mpiexec /home/VASP/vasp > "$PBS_JOBID".out.$j.$i rm WAVECAR CHG CHGCAR done done • for i in 250 350 450 500 • do • for j in 8 10 12 15 20 • do • cat >INCAR<<! • ENCUT = $i • ISMEAR = 0 • SIGMA =0.2 • NSW = 0 • EDIFF = 0.001 • ! • cat >POSCAR<<! • H2 • 1.0 • $j 0 0 • 0 $j 0 • 0 0 $j • 2 • Cart • 0 0 0 • 0 0 0.8 • ! • mpiexec /home/VASP/vasp > "$PBS_JOBID".out.$i.$j • rm WAVECAR CHG CHGCAR • done • done
How to move atoms? Energy is minimum. The gradient is zero.
Searching direction New start point Optimized step
Searching direction New start point New start point Optimized step Stop searching
设有二次函数:f(x) = 1/2 (x - x*)TA(x - x*) , 其中A是n×n对称正定矩阵,x*是一个定点,函数f(x)的等值面 1/2 (x - x*)TA(x - x*) = c 是以x*为中心的椭球面. 由于 ▽f(x*) = A(x - x*) = 0,A正定,因此x*是f(x)的极小点。设x(1)是在某个等值面上的一点,该等值面在点x(1)处的法向量▽f(x(1)) = A(x(1) - x*)。 设d(1)是这个等值面在x(1)处的一个切向量, d(1)与▽f(x(1))正交,即d(1)T▽f(x(1)) = 0; 另外,d(2) = x* - x(1), 有 d(1)TAd(2) = 0, 等值面上一点处的切向量与由这一点指向极小点的向量关于A共轭
Conjugate Gradient method Steepest descent method Conjugate Gradient
Structural relaxation To obtain the ground state relaxed geometry of the system. the equilibrium lattice constants a given ionic configuration the forces obtained these forces are greater than some minimum tolerance the ions are moved in the direction of the forces
IBRION-tag in VASP 0 Standard ab-initio molecular dynamics A Verlet algorithm is used to integrate Newton's equations of motion. POTIM supplies the timestep in femto seconds. 1 A quasi-Newton algorithm is used to relax the ions into their instantaneous ground state. 2 A conjugate-gradient algorithm i)In the first step ions :the direction of the steepest descent. The energy and the forces are recalculated. ii) Interpolation of the change of the total energy and of the forces, then a corrector step. iii) After the corrector step the forces and energy are recalculated and it is checked whether the forces contain a significant component parallel to the previous search direction. 1 initial position 2 trial step 3 corrector step, i.e. positions corresponding to anticipated minimum 4 trial step 5 corrector step …
Energies vs. Forces H2 molecular C -2.629 3.681 0.018 C -1.415 2.939 0.026 C -1.493 1.580 0.024 C -2.618 0.818 -0.011 C -3.851 1.541 0.016 C -3.904 2.981 -0.033 H -2.682 4.664 0.020 H -0.585 3.486 -0.047 H -0.544 1.065 -0.022 H -2.613 -0.192 -0.045 H -4.764 1.084 -0.040 H -4.683 3.483 0.032 C6H6 molecular
AWK and VASP aa=`awk '$2 == "TOTAL-FORCE" {print NR}' OUTCAR` a=`expr $aa + 1` bb=`expr $aa + $num + 2` awk 'NR > a && NR < bb {print $4,$5,$6}' a=$a bb=$bb OUTCAR >force.1 paste site.1 velocity.1 force.1 >all.1 awk '{print $1+$4*h+$7*h*h/2/m*0.0096,$2+$5*h+$8*h*h/2/m*0.0096, $3+$6*h+$9*h*h/2/m*0.0096}' h=$dt m=$mass all.1>site.2 paste velocity.1 force.1 force.2 >all.2 awk '{print $1+$4*h/2/m*0.0096+$7*h/2/m*0.0096, $2+$5*h/2/m*0.0096+$8*h/2/m*0.0096, $3+$6*h/2/m*0.0096+$9*h/2/m*0.0096}' h=$dt m=$mass all.2>velocity.2
Review of QN,SD and CG Conjugate Gradient Steepest descent method
Periodical structures three lines the three lattice vectors defining the unit cell of the system are given (first line corresponding to the first lattice vector, second to the second, and third to the third Latticeofgraphene & latticeofgrapheneD
K-mesh -- Monkhorst-Pack 与绝缘体比较,金属需要更密的 k-sampling to define the Fermi surface precisely。 q 取奇数,包括原点,不包括边界; q 取偶数,不包括原点,不包括边界; 。 Special points for Brillonin-zone integrations, PRB, 13,5188(1976)
KPOINTS test KPOINTS & IBZKPT • Automatic mesh • 0 ! number of k-points = 0 ->automatic generation scheme • Monkhorst-Pack ! select Monkhorst-Pack (first letter is significant) • Gamma • 5 5 1 ! size of mesh (4x4x4 points along b1, b2, b3) • 4 4 1 • 0. 0. 0. ! shift of the k-mesh
LREAL .FALSE. projection done in reciprocal space .TRUE. projection done in real space, This expression can be evaluated in reciprocal or real space: In reciprocal space the number of operations scales with the size of the basis set (i.e. number of plane-waves). In real space, the projection-operators are confined to spheres around each atom. Therefore the number of operations necessary to evaluate one does not increase with the system size (usually the number of grid points within the cut-off-sphere is between 500 and 2000). The number of atom <20 .FALSE. Otherwise, .TRUE.
Strain stabilized nanostructures PRL 100, 175503 (2008) ; 107, 236101 (2011)
Li on GaN surface mu=[-3 -2]; for ii=1:size(a,1) plot(mu,a(ii,2)-a(ii,1)*mu,'r') hold on end clear %Li a=[ 0 -.18783519E+03 1 -.76238692E+03/4 1.5 -.76726294E+03/4 2 -.77302221E+03/4 3 -.78305476E+03/4 %4 -.19702575E+03 4 -.19801252E+03 8 -.20551497E+03 ];
Defect formation energies NELECT @ INCAR Default (number of valence electrons) NELECT = number of electrons
Homework Review for the phase diagram. Sending to 17273799@qq.com when ready For lecture notes, refer to http://www.compphys.cn/~xbyang/