gromacs讀入碳納米管的方法

1、與gromacs讀入碳納米管的方法讀入碳納米管的方法有人問我leap怎么讀入碳納米管，這里說一下，順便也說一下gromacs讀入的方法。只說操作方法，不涉及選用力場等問題。碳納米管結(jié)構(gòu)比較大，而結(jié)構(gòu)很有規(guī)律，一般僅有C原子。若不考慮邊界問題，可以將全部電荷當(dāng)成0，這種近似并不會對結(jié)果真實性有多少影響。對于碳納米管這樣的特征，適合當(dāng)成一般的小分子來處理，但不宜用antechamber而且也沒必要用。先用Nanotubemodeler生成一

2、段碳納米管，含240個C原子，導(dǎo)出結(jié)構(gòu)為tube.pdb。但是會發(fā)現(xiàn)，原子名都是相同的，而leap要求結(jié)構(gòu)中相同殘基中每個原子都有獨立的名字，所以需要先改名。比較方便的方法是用ultraedit，打開后開啟列模式，選定從第1個原子到最后一個原子的15和16列(即原子名C后面的兩列)，選ColumnNumber，OK，即把原子名改為了C1至C240，保存為tubename.pdb。這里假定C原子使用gaff力場的c2原子類型對應(yīng)的力場參數(shù)

3、，依次運行xleapfleaprc.gaffa=loadpdbsobtubename.pdbbondbydistancea默認(rèn)是2埃內(nèi)成鍵，故碳納米管的碳原子都會正常成鍵，可用edita檢查。然后把碳納米管內(nèi)全部240個原子都設(shè)成c2原子類型，可以用批處理腳本。由于leap只能一條一條運行批處理文件中的內(nèi)容而不支持shell腳本的循環(huán)，需要將循環(huán)轉(zhuǎn)化成普通命令腳本再在leap運行。寫腳本：f((i=1ileapdodone在shell下

4、運行，得到leapdo，然后在leap里運行sourceleapdo即可checka進(jìn)行參數(shù)檢查，應(yīng)該OK。如果想改力場參數(shù)或缺少某些參數(shù)，可以editgaff往里添加和修改。如果要加溶劑和普通情況無異，這里略過。最后saveamberparmatube.toptube.inpcrd=============gromacs讀入碳納米管還是用上面的tube.pdb為例子。小分子在gromacs中常用prodrg來處理，但是理由和antec

5、hamber一樣并不適合，而且有更方便的辦法。gromacs提供了一個構(gòu)建拓?fù)湮募墓ぞ選2top，專適用于構(gòu)建結(jié)構(gòu)規(guī)律性很強的體系。x2top中有一些bug和“規(guī)矩”，而且在不同gromacs版本中bug和“規(guī)矩”還不一樣，這里使用gromacs4.0.4的x2top。首先需要寫n2t文件，比如使用ffgmx2力場，就在力場文件所在文件夾里寫一個ffgmx2.n2t，內(nèi)容如下：CCX0.012.0113C0.14C0.14C0.14C

6、CX0.012.0112C0.14C0.14CCX0.012.0111C0.14第一行說明如果體系中任何一個C原子，與周圍3個C原子的距離都在0.14nm左右（判斷標(biāo)準(zhǔn)大概為正負(fù)10%），就把它當(dāng)作CX原子類型（原子類型名字自定，可以不屬于力場包含的原子類型），電荷為0.0，相對原子質(zhì)量為12.011。第二行和第三行與第一行意義類似，即代表與兩個C原子和與一個C原子成鍵的C原子被當(dāng)成什么類型、多少電荷、多少原子質(zhì)量，此例中將碳納米管所有

7、原子都當(dāng)成CX。運行x2topftube.pdbotube.topffgmx2這樣就得到tube.top，里面包含這個碳納米管全部鍵結(jié)項。ff代表用什么力場，此處即ffgmx2，如果ff則是出現(xiàn)列表自行選擇。這里如果使用ffgmx力場即便正確寫了ffgmx.n2t也可能無法正確執(zhí)行，即不能正確根據(jù)距離判斷成鍵，此時應(yīng)換用別的力場來執(zhí)行x2top。如果運行后程序卡住不動，可嘗試加nopbc解決，這是一個bug。默認(rèn)情況下，x2top會自動

8、在.top里相應(yīng)鍵結(jié)項加上力場參數(shù)，這種方式自動加入的力場參數(shù)并不是根據(jù)力場中相應(yīng)原子類型間的鍵結(jié)參數(shù)加入的。平衡距離鍵角二面角就是tube.pdb中的相應(yīng)項的距離鍵角二面角，鍵的力常數(shù)鍵角力常數(shù)二面角力常數(shù)默認(rèn)是4000004005，可以分別用kb、kt、kp設(shè)定。如果不想讓其自動加入，可以加上noparam。如果n2t中把C原子轉(zhuǎn)換為ffgmx2力場中已有的類型，比如CB，那么此時這個.top搭配tube.pdb已經(jīng)可以直接用于模擬

9、了。此例將C設(shè)為力場中沒有的原子類型CX是為了假設(shè)需要加入專門的力場參數(shù)的情況，此時需要定義新的原子類型。首先加入VDW和14VDW參數(shù)，分別在ffgmx2nb.itp里的[atomtypes]和[nonbond_params]里面添加。如果前面x2top時加了noparam，還需要設(shè)定鍵結(jié)參數(shù)。分兩種情況：1直接在.top里面每一個鍵結(jié)項后面加入相應(yīng)力場參數(shù)，類似于x2top自動在.top里添加力場參數(shù)的方式，用ultraedit等軟