计算化学公社

 找回密码 Forget password
 注册 Register
Views: 2705|回复 Reply: 5

[Quantum ESPRESSO] QE中批量固定原子的小脚本QEfix.sh

[复制链接 Copy URL]

327

帖子

9

威望

1995

eV
积分
2502

Level 5 (御坂)

发表于 Post on 2022-2-8 16:31:03 | 显示全部楼层 Show all |阅读模式 Reading model
本帖最后由 丁越 于 2022-8-6 00:19 编辑

QE中批量固定原子的小脚本QEfix.sh


      虽然在pwgui中提供了固定原子的选项,但是一个一个地固定原子坐标显得颇为麻烦,遂写了个小脚本可以很方便的固定原子坐标。脚本的使用很简单,例如将rutile的结构载入gview中,然后选择所要固定的原子后,在tools -> atom selection中就可以复制所要固定原子的序号了。接下来,通过pwgui生成如结构优化输入文件rutile.relax.in(不需要固定原子),然后执行 ./QEfix.sh rutile.relax.in 后,脚本首先会提示需要固定原子的序号,粘贴回车;紧接着会提示如何固定,如输入 0 0 0就表示将x、y、z三个方向同时固定。
      接下来说说脚本的思路:原子的序号样式是这样的 2,3,4-10,因此,首先我们需要将这些序号拆解成一个一个的单个序号形式(2,3,4,5....10),然后写个for循环在批量去在每个需要固定的原子后面添加如0 0 0。拆解原子序号时,我们先用awk以逗号为分割符将其拆开,然后存进第一个临时文件tmp1.dat;接下来,我们需要判断tmp1.dat中数据假如是单个数字的,那就直接把它存进临时文件tmp2.dat中,假如是以4-10这样的数据,我们设置Min=4,Max=10,然后利用for循环把Min~Max中的数据缀加到tmp2中,此时,整个原子序号就被拆解成了单个序号存进了tmp2.dat中。接下来,由于输入文件中第一个原子所在行数不是1,那么我们还需要将tmp2.dat数据都加上行数得到tmp3.dat; 最后,再次利用for循环按照tmp3.dat中原子所在行数去添加如0 0 0。脚本中,[ "$i" =~ ^[0-9]+$ ]是比较常用的判断字符串包含关系的语法。sed -i "${i}s/$/${if_pos}/g" $1中$表示行尾,s表示替换,g为global,$1意思是脚本文件接受的第一个参数。

  1. #!/bin/bash
  2. #usage: QEfix.sh [qe input file]

  3. fname1=`echo $1`
  4. echo ' Input indices of the atoms to be constraint (fixed), e.g. 1,5,9-12,14-18: '
  5. read num
  6. echo ' Input if_pos(i), e.g. 0 0 0, which means atomic coordinates at x, y, and z direction will be fixed simutaneously:'
  7. read if_pos
  8. echo

  9. pos=`grep -n ATOMIC_POSITIONS ${fname1} |cut -d ":" -f 1 |head -1`
  10. echo $num|awk -F , '{for(i=1;i<=NF;i++) print $i}' > ${prefix}_tmp1.dat

  11. #Generating the atoms Index into tmp2.dat
  12. for i in $(cat ${prefix}_tmp1.dat);do
  13.         if [[ "$i" =~ ^[0-9]+$ ]];then
  14.                 echo $i >> ${prefix}_tmp2.dat
  15.         else
  16.                 Min=$(echo $i |awk -F \- '{print $1}')
  17.                 Max=$(echo $i |awk -F \- '{print $2}')
  18.                 for j in $(seq $Min $Max)
  19.                 do
  20.                        echo $j >> ${prefix}_tmp2.dat
  21.               done               
  22.         fi
  23. done

  24. #fixing atoms in the QE input file
  25. for i in $(cat ${prefix}_tmp2.dat)
  26. do
  27.         line=$(($i + $pos))
  28.         echo $line >> ${prefix}_tmp3.dat
  29. done

  30. for i in $(cat ${prefix}_tmp3.dat)
  31. do
  32.         sed -i "${i}s/$/    ${if_pos}/g" ${fname1}
  33. done

  34. rm -f ${prefix}_tmp*


复制代码






QEfix.sh

1021 Bytes, 下载次数 Times of downloads: 5

评分 Rate

参与人数
Participants 1
威望 +1 收起 理由
Reason
sobereva + 1

查看全部评分 View all ratings

自由发挥,野蛮生长

3015

帖子

3

威望

1万

eV
积分
14751

Level 6 (一方通行)

从头算界孔乙己

发表于 Post on 2022-2-8 19:30:12 | 显示全部楼层 Show all
本帖最后由 卡开发发 于 2022-2-8 19:51 编辑

其实利用ase gui可以直接在ase自带的图形界面直接批量选择原子进行固定:
(1)先启动ase
  1. $ ase gui CeO2.xsd
复制代码
(2)然后界面上工具栏的Tools->Constraints调出固定的菜单。
(3)此时可以通过鼠标框选原子(被选中的原子的轮廓会加粗).
(4)再到Constraints呼出的菜单点击Fix按钮,此时选中的原子上会显示X记号表面这些原子被固定了。
(5)保存,文件名设置为*.pwi就行,ase能识别为qe的in文件(小心不要用poscar之类的命名,否则会被误认为vasp的结构)。

其实磁矩也可以使用ase的gui设置,大家可以自己试试。最后,即便导出的结果不能满意,其实也可以按照ase的json格式存储,再简单写一些程序进一步处理。事实上在上述基础上之前我用PyQT设计出一版额外增加QE参数的图形界面,目前不方便公开,因为一方面是和别人合作的项目,另一方面就是QE的输入方式在工程角度看稍有一些缺陷,很多东西也不是太完善。


评分 Rate

参与人数
Participants 1
eV +3 收起 理由
Reason
hebrewsnabla + 3

查看全部评分 View all ratings

近期不及时回复。欢迎无偏见非商业的学术讨论,但是看家本领和课题组的传统艺能别人会毫无保留告诉你?

327

帖子

9

威望

1995

eV
积分
2502

Level 5 (御坂)

 楼主 Author| 发表于 Post on 2022-2-8 20:24:26 | 显示全部楼层 Show all
卡开发发 发表于 2022-2-8 19:30
其实利用ase gui可以直接在ase自带的图形界面直接批量选择原子进行固定:
(1)先启动ase(2)然后界面上 ...

谢谢卡开发发老师,ASE我还没用过,主要是不会python过段时间再学下
自由发挥,野蛮生长

3015

帖子

3

威望

1万

eV
积分
14751

Level 6 (一方通行)

从头算界孔乙己

发表于 Post on 2022-2-8 20:41:14 | 显示全部楼层 Show all
丁越 发表于 2022-2-8 20:24
谢谢卡开发发老师,ASE我还没用过,主要是不会python过段时间再学下

python可以学学,一些轻量化的处理会轻松很多,ase也就是其中一个生态资源而已。
近期不及时回复。欢迎无偏见非商业的学术讨论,但是看家本领和课题组的传统艺能别人会毫无保留告诉你?

100

帖子

0

威望

3380

eV
积分
3480

Level 5 (御坂)

发表于 Post on 2022-8-5 09:43:02 | 显示全部楼层 Show all
感谢您的脚本!似乎有个微小问题。固定参数(例如1 1 1)与原子坐标z值相连,没有空格。
  1. N     6.690430   15.690897    3.6265201 1 1
复制代码

可以将脚本第36行作如下调整
  1.        #original
  2.         sed -i "${i}s/$/${if_pos}/g" $1
  3.        #modified
  4.         sed -i "${i}s/$/ ${if_pos}/g" $1
复制代码

即可。
gcq#413: "I know poetry is not dead, nor genius lost; nor has Mammon gained power over either, to bind or slay; they will both assert their existence, their presence, their liberty and strength again one day." (Jane Eyre in Jane Eyre by Charlotte Bronte)

327

帖子

9

威望

1995

eV
积分
2502

Level 5 (御坂)

 楼主 Author| 发表于 Post on 2022-8-6 00:11:34 | 显示全部楼层 Show all
qczgzly 发表于 2022-8-5 09:43
感谢您的脚本!似乎有个微小问题。固定参数(例如1 1 1)与原子坐标z值相连,没有空格。

可以将脚本第36 ...

谢谢提醒,我刚刚再写QEtoolkit2.0版本的时候也发现了这个问题,估计明后天就发出来了,用它生成QE的输入文件将会非常方便。
自由发挥,野蛮生长

本版积分规则 Credits rule

手机版 Mobile version|北京科音自然科学研究中心 Beijing Kein Research Center for Natural Sciences|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949号-1 )|网站地图

GMT+8, 2023-2-7 02:31 , Processed in 0.391468 second(s), 26 queries .

快速回复 返回顶部 返回列表 Return to list