仓库源文，站点原文

layout: post title: 统计轨迹中分子速度大小沿某一方向的分布 categories:

科 tags:
awk
gmx

2017年03月29日 09:25:09

写过一个简单的awk脚本, 用来统计轨迹中水分子的质心速度大小沿坐标的变化情况, 目的是可以用来拟合计算粘度. 具体的计算方法可参考:

A Comparison Of The Value Of Viscosity For Several Water Models Using Poiseuille Flow In A Nano-channel

A. P. Markesteijn, Remco Hartkamp, S. Luding, J. Westerweel; J. Chem. Phys. 136(13):134104, 2012; 10.1063/1.3697977

脚本只能处理gro格式的轨迹, 效率不是很高, 但基本能用. 放在这里, 存个档, 也供需要的人参考吧.

<table class="highlighttable"><th colspan="2" style="text-align:left">zvx.bsh</th><tr><td><div class="linenodiv" style="background-color: #f0f0f0; padding-right: 10px"><pre style="line-height: 125%"> 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24</pre></div></td><td class="code"><div class="highlight" style="background: #f8f8f8"><pre style="line-height: 125%">awk ' BEGIN { mo=16; mh=1; M=mo+2mh # 定义原子质量 Ntot=10 # 定义要计算的帧数 Zmin=3; Zmax=26; dZ=0.2 # 定义计算区间, 分格间距 Nfrm=0; N=int((Zmax-Zmin)/dZ) for(i=0; i<=N; i++) { Z[i]=Zmin+idZ; Vx[i]=0; Nwat[i]=0 } }

NF==1 { Nfrm++ } NF>6 && $2~/OW/ { # 根据 OW 原子统计 z = mo$(NF-3); vx = mo$(NF-2) # O 的z坐标, x方向速度 getline; z += mh$(NF-3); vx += mh$(NF-2) # H1的 getline; z += mh$(NF-3); vx += mh$(NF-2) # H2的

z /<span style="color: #666666">=</span> M; vx /<span style="color: #666666">=</span> M    <span style="color: #008800; font-style: italic"># 质心坐标和速度</span>
<span style="color: #B8860B">i</span><span style="color: #666666">=</span>int<span style="color: #666666">((</span>z-Zmin<span style="color: #666666">)</span>/dZ<span style="color: #666666">)</span>
Nwat<span style="color: #666666">[</span>i<span style="color: #666666">]</span>++; Vx<span style="color: #666666">[</span>i<span style="color: #666666">]</span> +<span style="color: #666666">=</span> vx
<span style="color: #AA22FF; font-weight: bold">if</span><span style="color: #666666">(</span>Nfrm&gt;Ntot<span style="color: #666666">)</span> <span style="color: #AA22FF">exit</span>

} END { print "# Z: ["Zmin":"Zmax":"dZ"] Frames: " Nfrm for(i=0; i<=N; i++) if(Nwat[i]>0) printf"%f %f\n", Z[i], Vx[i]/Nwat[i] } ' File.gro >Z-Vx.xvg </pre></div> </td></tr></table>

2017-04-19 21:30:16 人在囧图 李老师：你好请问下，我跑完gromacs后得到的速度和坐标信息在.trr文件中，您这个脚本只能识别.gro格式的文件，怎么可以有效的把.trr格式的文件转化为可识别的.gro文件？谢谢李老师。
2017-04-21 04:28:59 Jerkwin 用 trjconv -f -n -o File.gro 就可以得到gro文件了, 注意其中水的氧原子名称必须是 OW, 否则的话, 你要修改脚本

仓库源文，站点原文

评论