將GROMACS的原子電荷信息讀入VMD的方法
2023-Jun-8注：推薦用http://bbs.keinsci.com/thread-37839-1-1.html里介紹的VMD載入tpr文件的插件，比本文的做法更方便。VMD載入tpr文件后直接就有了tpr中的原子電荷信息。
將GROMACS的原子電荷信息讀入VMD的方法
文/Sobereva @北京科音
First release: 2017-Mar-27 Last update: 2022-Apr-12
《使VMD實時顯示gromacs軌跡中原子的受力》（http://www.shanxitv.org/36）
《使VMD讀入Gromacs產生的trr軌跡中速度信息的方法》（http://www.shanxitv.org/117）
有很多VMD里面的插件，比如顯示偶極矩、計算靜電勢、繪制紅外光譜的插件都依賴于原子電荷，有時候我們還需要編寫依賴于原子電荷的腳本、使用依賴于原子電荷的選擇語句，因此把GROMACS模擬時用的原子電荷載入VMD也是非常重要的，這里說說怎么做。
首先在這里下載筆者開發的gmxoutchg程序：http://www.shanxitv.org/soft/gmxoutchg_1.1.rar
其中gmxoutchg.exe是Windows版可執行文件，沒后綴的gmxoutchg是Linux下可執行文件，gmxoutchg.f90是源文件。經測試此程序兼容GROMACS 2016.1、2018.8和2019.3版，對其它版本兼容性未經測試。
tpr文件里蘊含了模擬所需的一切信息，包括原子電荷。用GROMACS里的dump命令將二進制的tpr文件轉化成文本文件dump.txt，執行以下命令：
gmx dump -s test.tpr > dump.txt
然后將dump.txt放到gmxoutchg.exe所在目錄下，運行gmxoutchg.exe，程序就會解析其中的數據，在當前目錄下產生charge.txt。其中包含體系所有原子的原子電荷，順序和結構文件里的原子順序完全一致。其中第一列、第二列、第三列分別是原子電荷、分子序號、原子序號。
將charge.txt放到VMD目錄后，在載入對應的結構/軌跡后，可以用以下tcl腳本將原子電荷從charge.txt中讀入。
set sel [atomselect top all]
set natom [$sel num]
set rdchg [open "charge.txt" r]
set chglist {}
for {set iatm 0} {$iatm<=[expr $natom-1]} {incr iatm} {
gets $rdchg line
scan $line "%f" chg
lappend chglist $chg
}
$sel set charge $chglist
close $rdchg
如果想驗證是否正確載入了，可以用[atomselect top all] get charge命令把所有原子的原子電荷輸出出來（當然，也可以把all換成選擇語句只輸出指定部分的原子電荷）。
也可以將Coloring method設成Charge，直接從顏色上檢驗是否正確載入了。下面是乙醇盒子，在默認色彩刻度(RWB)下越紅代表原子電荷越負，越藍代表原子電荷越正。
順帶一提，對于一些體系通過恰當設定顯示方式，可以令原子電荷分布一目了然。比如下圖的分子，一個顯示方式設為Licorice并且把鍵調細，另一個顯示方式是CPK，把圓球調大，用透明材質，用charge來著色，效果挺不錯。
有了原子電荷信息可以直接使用VMD的依賴于原子電荷的插件了，比如Extensions-Visualization-Dipole Moment Watcher可以觀看基于原子電荷計算的偶極矩矢量，就是上圖的紅色箭頭。
PS：還有一種把原子電荷載入的方法是使用
gmx editconf -f nico.tpr -mead nico.pqr
這樣nico.tpr里對應的原子電荷就被輸出到了nico.pqr的電荷那一列了，pqr文件里還有原子半徑信息，把nico.pqr載入VMD時也會讀入其中的原子電荷。但這種做法致命缺點是原子電荷保留精度太低，只有兩位小數，當涉及一些靜電作用的分析時誤差可能較大。