| 查看: 4482 | 回复: 8 | |||
valenhou001至尊木虫 (职业作家)
|
[交流]
从VASP的OUTCAR提取xyz以及结构优化的轨迹 已有4人参与
|
|
将以下内容考虑至一个文件，比如命名为getvaspxyz.sh，然后将它的属性改为可执行的, chmod +x getvaspxyz.sh。另外在linux下，还可能需要使用命令dos2unix来去掉文件里每行多余隐藏的字符：dos2unix getvaspxyz.sh。在含有OUTCAR文件的目录下运行getvaspxyz.sh，就可以取出原子的位置，以及原子所受的力。按xyz的格式输出。 xyz格式的文件可以用Jmol等软件来可视化。 一个缺点就是，如果所计算的体系所包含的原子个数比较多，另外结构弛豫的步数比较多的话，在得到轨迹文件比较慢。 #!/bin/bash pwd1=`echo ${PWD##*/}` natom=`grep \'NIONS\' OUTCAR|tail -1 |awk \'{printf "%s\\n", $12}\'` last=`sed -n \'/total drift:/=\' OUTCAR |tail -1` iaf=`echo "$last - 2 - $natom" |bc -l ` echo $natom > opt.xyz sed -n "$iaf",+"$natom"p OUTCAR >>opt.xyz echo "The latest atomic positions are written to opt.xyz!" iptc=`grep \'ions per type =\' OUTCAR|tail -1 ` ipt=(`echo "${iptc##*=}"`) ats=(`grep \'VRHFIN =\' OUTCAR | sed -e \'s/=/ /g\' |sed -e \'s/:/ /g\'| awk \'{printf "%s\\n", $2}\' `) ntp=${#ipt[@]} i=2 for it in `seq 0 $((${ntp}-1))`;do for j in `seq 1 ${ipt[$it]}`;do i=`echo "$i +1 " |bc` elem=${ats[$it]} sed -i "$i s/^/$elem/" opt.xyz done done cp opt.xyz opt-"$pwd1".xyz nstep=(`sed -n \'/total drift:/=\' OUTCAR `) etot=(`grep \'energy without entropy=\' OUTCAR |awk \'{printf "%20.8f\\n", $4}\' `) myFile="trajec.xyz" if [ ! -e "$myFile" ]; then echo $natom > $myFile sed -i \'1d\' $myFile fi for k in `seq 0 $((${#nstep[@]}-1))`;do iaf=`echo "${nstep[k]} - 2 - $natom" |bc -l ` if [ $k -eq 0 ];then echo $natom > $myFile else echo $natom >> $myFile fi sed -n "$iaf",+"$natom"p OUTCAR >> $myFile i=`echo "$k * ($natom+2) + 2"|bc -l` sed -i "$i s/^/#Etotal: ${etot[k]} eV/" $myFile for it in `seq 0 $((${ntp}-1))`;do for j in `seq 1 ${ipt[$it]}`;do i=`echo "$i +1 " |bc` elem=${ats[$it]} sed -i "$i s/^/$elem/" $myFile done done done echo "The atomic positions in each optimization step are written to trajec.xyz!" cp trajec.xyz trajec-"$pwd1".xyz |
回复此楼» 本帖已获得的红花（最新10朵）
lion_tiger» 猜你喜欢
为什么国自然不能直接公布
已经有5人回复
-
梦想
已经有7人回复
-
学科评审组评审是指会评吗？
已经有4人回复
-
面上合作单位盖章
已经有8人回复
-
面上函评意见出来了，像什么等级？
已经有15人回复
-
为什么资助数各大高校都创新高，自己申请怎么就这么难
已经有14人回复
-
小白求助 投论文要求的highlights应该如何写
已经有4人回复
-
麻烦专家们看看评委们的意见（F口面上）
已经有12人回复
-
怎么看青基中了没有啊
已经有6人回复
-
国社科又开始会评了，不知道这次命运如何
已经有11人回复
lion_tiger
至尊木虫 (正式写手)
- 应助: 23 (小学生)
- 金币: 16647
- 红花: 9
- 帖子: 579
- 在线: 1525.3小时
- 虫号: 1598434
- 注册: 2012-02-03
- 性别: GG
- 专业: 固体力学
2楼2019-10-29 16:45:23
zuocuiping
木虫 (职业作家)
- 应助: 14 (小学生)
- 金币: 4975.9
- 散金: 966
- 红花: 60
- 帖子: 4045
- 在线: 1600.3小时
- 虫号: 890422
- 注册: 2009-11-01
- 性别: GG
- 专业: 半导体光电子器件
3楼2019-10-29 18:01:05
valenhou001
至尊木虫 (职业作家)
- 1ST强帖: 13
- 应助: 241 (大学生)
- 金币: 25701.7
- 散金: 602
- 红花: 166
- 帖子: 3782
- 在线: 873.8小时
- 虫号: 1007127
- 注册: 2010-04-27
- 专业: 凝聚态物性 II ：电子结构
4楼2019-10-29 20:13:04
zuocuiping
木虫 (职业作家)
- 应助: 14 (小学生)
- 金币: 4975.9
- 散金: 966
- 红花: 60
- 帖子: 4045
- 在线: 1600.3小时
- 虫号: 890422
- 注册: 2009-11-01
- 性别: GG
- 专业: 半导体光电子器件
5楼2019-10-29 22:01:55
hy1987007
木虫 (著名写手)
- 应助: 12 (小学生)
- 金币: 20203.7
- 散金: 2
- 红花: 2
- 帖子: 2555
- 在线: 459小时
- 虫号: 1560286
- 注册: 2012-01-02
- 性别: GG
- 专业: 凝聚态物性 II ：电子结构
6楼2019-10-29 23:30:30
轻轻飘过之
金虫 (正式写手)
- 应助: 5 (幼儿园)
- 金币: 2897.8
- 散金: 189
- 红花: 2
- 帖子: 474
- 在线: 217.3小时
- 虫号: 1401636
- 注册: 2011-09-14
- 性别: MM
- 专业: 凝聚态物性 II ：电子结构
7楼2019-10-30 09:35:05
valenhou001
至尊木虫 (职业作家)
- 1ST强帖: 13
- 应助: 241 (大学生)
- 金币: 25701.7
- 散金: 602
- 红花: 166
- 帖子: 3782
- 在线: 873.8小时
- 虫号: 1007127
- 注册: 2010-04-27
- 专业: 凝聚态物性 II ：电子结构
8楼2019-11-13 10:10:24
valenhou001
至尊木虫 (职业作家)
- 1ST强帖: 13
- 应助: 241 (大学生)
- 金币: 25701.7
- 散金: 602
- 红花: 166
- 帖子: 3782
- 在线: 873.8小时
- 虫号: 1007127
- 注册: 2010-04-27
- 专业: 凝聚态物性 II ：电子结构
|
在用jmol显示这个脚本所提取出来的轨迹时。通过在jmol的控制台使用如下命令： set vectorScale 5.0 set vectorsCentered set vectorSymmetry FALSE set vibrationPeriod 2.0 set vibrationScale 2.0 set defaults PyMOL vector ON vector 0.05 vibration ON 可以用来显示原子所受的力的矢量。 或者将以上几行命令放到一个文件比如名为"label.spt"的文件中，然后将该文件通过鼠标拖到jmol已打开了原子轨迹显示的窗口，那这些原子上就会显示出相应所受的力的矢量，以及原子沿着这个矢量上的振动，这样可以很清楚地看到那个原子受的力是很大。 |
9楼2019-11-13 10:15:23
20