關於towhee 回答zhang_jaj的問題-轉載小木蟲

如下是zhang_jaj 的問題:
你好，我想用towhee作NVT下argon氣在surface面上nucleation過程，但願你給一些指導，首先我想把potential選爲LJ，surface與particles也是LJ做用，
1）問題是這個surface是須要用particles來make嗎，若是是的話那麼構成它的particles該怎麼固定
2 ）或者，直接用towhee裏面的LJ9-3 wall嗎，想知道這種wall是與particles怎麼做用的，它不是用小粒子構成的嗎？
3）模擬完成後，液體與固體的結構判斷該怎麼作，有什麼程序代碼嗎？
4）這個nucleation的rate該怎麼算（是單位時間單位體積內產生臨街晶核的個數），我想知道程序中怎麼實現，我沒發現有那個參數相關的，或者是須要本身來加東西嗎？
5）還有這個nucleation barrier該怎麼來實現？
6）towhee不能並行，那麼其實一個任務就是在一個process上處理的；若是個人我的電腦比服務器的主頻高，那麼是否是用臺式機比服務器更快呢（由於臺式機是4核的，它是否是能本身把towhee的任務來用4核一塊兒處理？）？
7）但願您能給點關於towhee代碼學習的經驗，就是對於我要作的這個方向的。

這是個人回答:
實在不敢當, 我如今也是學習階段, 由於沒作過相似的體系, 對於你的問題只能有一些建議.
1）問題是這個surface是須要用particles來make嗎，若是是的話那麼構成它的particles該怎麼固定
可使用particle模型來模擬, 固定比較簡單, 在move設置中不讓這些粒子移動就好了.
2 ）或者，直接用towhee裏面的LJ9-3 wall嗎，想知道這種wall是與particles怎麼做用的，它不是用小粒子構成的嗎？
這個我之前沒有用過, 不過文獻當中確實有人這樣作. 這種勢能的使用狀況能夠參考例子Fris_Walls, 另外上網搜了一下, 這幾個帖子和文章可能對你有用.
http://jcp.aip.org/resource/1/jcpsa6/v67/i5/p2384_s1
http://jcp.aip.org/resource/1/jcpsa6/v73/i8/p4050_s1
http://lammps.sandia.gov/threads/msg24356.html
http://www.sklogwiki.org/SklogWi ... ard-Jones_potential
這種勢是粒子勢能轉化而來的, 第四個參考文獻有推導過程, 主要的優點是能夠大大提升速度. lammps 的帖子跟你的狀況有點像, 使用MD作的, 能夠借鑑.
3）模擬完成後，液體與固體的結構判斷該怎麼作，有什麼程序代碼嗎？
這個主要是分析RDF吧, 代碼在 utils 中有, 能夠根據你的狀況進行修改.
4）這個nucleation的rate該怎麼算（是單位時間單位體積內產生臨街晶核的個數），我想知道程序中怎麼實現，我沒發現有那個參數相關的，或者是須要本身來加東西嗎？
這個是須要後處理才能獲得的吧. towhee自己不是計算成核過程的. 應該須要編寫代碼實現. 不過MC 如何計算時間呢? 若是用MC步數是否是會有問題? 是否是用MD更合適一點?
5）還有這個nucleation barrier該怎麼來實現？
這個應該是進行不一樣溫度點的模擬, 做速率-溫度關係圖, 而後經過擬合熱力學公式獲得吧.
6）towhee不能並行，那麼其實一個任務就是在一個process上處理的；若是個人我的電腦比服務器的主頻高，那麼是否是用臺式機比服務器更快呢（由於臺式機是4核的，它是否是能本身把towhee的任務來用4核一塊兒處理？）？
towhee的不能用多個核並行計算一個任務, 而事實上也不須要. 通常來講, 模擬多個任務能夠等價於模擬一個大任務, MC的本質決定它是這樣的性質, MD可不行. 因此儘可能仍是使用towhee提供的並行模式. 注意選擇不一樣的隨機種子.
7）但願您能給點關於towhee代碼學習的經驗，就是對於我要作的這個方向的。
關於代碼, 其實我也沒有太多經驗. 我主要是對結果進行後處理. towhee 自己是一個general Monte Carlo code, 不是專門爲某些體系設計的, 可是對於你的體系應該是能夠勝任的. 主要針對體系的代碼, 差異也主要是建模和分析結果兩個過程, 核心代碼都是相似的. 我我的認爲若是時間有限, 不要陷到源代碼中去.固然可以多瞭解一些代碼對於軟件的使用也是很是有好處的.
towhee 的代碼因爲須要我簡單的讀過一些, 仍是比較清楚的. 能夠嘗試從整體上去學習下.
對於結果的後處理的代碼, 建議從修改towhee 提供的代碼開始. 使用語言以本身最熟悉的爲準, 若是都湊合, 推薦 fortran90.
多學習 linux bash 腳本也是頗有好處的, 我如今簡單的數據處理基本都用它.
最後強調, 不要陷到代碼中去, 可以完成本身的用途就能夠了, 越簡單越直接越好. 多多思考跟你體系相關的問題, 多查文獻.

///但願能多多向你學習towhee知識，導師無論，本身一我的摸索。擔憂畢業十大問題啊，謝謝了
你們基本都差很少, 堅決信心, 能夠作的很好. 其實模擬計算這個方向比起作實驗仍是比較透明的, 重複文獻通常沒有問題, 技術問題通常論壇上都會回答.
我也是在學習的過程, 你們多多討論, 相互學習.

————————————————————————————————————————————————————————————————————

你好，想再請教下您。關於towhee參數設置問題：
1）.pm*的值是否是按照從小到大進行排列的呢。example中有個是
pmtraat 1.0 

pmtracm 0.67  

pmrotate 1.0
這種狀況，當要判斷是進行哪一種move時，是否是隻進行pmtraat呢？

2）.給體系加上一個LJ-9-3 wall時，裏面有一個參數是ljfrho（the number density of atoms in the integrateed wall potential）是wall上與體系內分子做用的點密度嗎？還有ljfsig與ljfeps都是人爲的指定大小嗎？

很是感謝！

————————————————————————————————

1. 是的，嘗試機率逐漸遞增。例如
pm1 0.1
pm2 0.5
pm3 1.0
pm4 1.0
那麼pm1～4的機率分別爲0.1，0.4，0.5，0. 只要出現1.0，那麼後面的pm機率就都爲0了
2.不是，參數的具體意義和設置看我發給你的幾個連接。

————————————————————————————————————————

你好，還想向您請教個問題關於hmatrix參數的。mannul上說這個座標矩陣能夠不正交，可是必須是right hand rule的。我在想若是用其餘非正交的座標矩陣時，程序在運行過程當中是否是速度就低了不少？？？（由於我的認爲正交座標裏各個原子間進行距離計算時要方便的多（如我在一個面上固定一個單層石墨烯））
還有個問題是，對於NVT系宗，若是我先讓N=1000，T=50K後，該怎麼選取box的大小呢（氣體），這時不一樣的box大小應該就對應不一樣的P吧。是利用PV=NRT嗎？這個box的大小是依據我想要的P用那個公式來調節嗎？，仍是normal atmosphere下算個大約的box大小？

很是感謝！！！

————————————————————————————————————————————————

沒有用過非正交的盒子，不過這個對速度影響應該不大吧。你能夠試試。
對於氣體的話可使用PV=NRT來估算，注意單位。壓力是根據你的體系要求指定的吧。

————————————————————————————————————————

還有關於LJ-9-3wall裏面的ljfrho參數，不知道你涉及沒有，這個我實在是理解不了是怎麼算的，用哪些原子，除以多大的空間（和做用範圍有關嗎）？由於英語很菜，對那句話理解不清楚，很是感謝！

————————————————————————————————————————

我看了一下LJ93wall的參數手冊，lifrho是指每立方埃內的分子數。好比你的wall是Si，每一個晶胞中分子數除以晶胞體積。若是你的wall是SiO2，可能就要創建位置重疊的兩個wall，Si和O分別創建wall，這樣的數密度分別是SiO2晶胞中Si的數目除以晶胞體積，O的數目除以晶胞體積。

——————————————————————————————————————————————

這個看了仍是不太明白。按照mannul應該是把wall當作一種通常性的surface，而這個ljfrho就是這個面上的點分佈密度。還有下面的那個sigma和epislon都是這個wall與體系中粒子的做用參數。（這個參數該怎麼定義呢），還有若是把這個wall當作Si之類的來計算ljfrho，那麼這些Si原子的排列該怎麼肯定(或者說這是一個無序的surface嗎)？你以前給的那個文章連接我一直下不了，不知道你有嗎，能夠給我一份嗎（zhjjaj@163.com）？關於LJ9-3的文章很是感謝！http://jcp.aip.org/resource/1/jc ... _s1?isAuthorized=no

————————————————————————————————————————————————

我從頭講好了。
咱們通常模擬時使用的是全粒子勢能好比UFF，amber等等，這種勢能每一個原子做爲獨立的粒子，有本身的力場參數。這樣在模擬時須要計算每一對粒子的相互做用能。若是咱們的體系是氣體或者液體分子和一個表面的相互做用，固然也可使用全粒子勢能參數去模擬，可是咱們也能夠把表面全部粒子的力場參數積分使其變成一個總體。9-3 wall 勢能是從 12-6 粒子勢能進行空間積分獲得的。積分推導的過程見：http://www.sklogwiki.org/SklogWi ... ard-Jones_potential , 建議本身能夠簡單推導一下。這樣作的好處是大大減小了粒子的數目（表面當成了一個具備特定形狀的粒子），從而減小了計算量。推導過程當中，原來12-6中的參數sigma ，epsilon意義沒有變化，直接繼承下來，另外新加入了一個參數ljfrhho，和體系的組成相關，即體相單位體積內的粒子數目。
我只是以Si和SiO2 爲例，你要根據你的體系去獲得sigma epsilon和ljfrhho。請問你的體系中wall的組成是什麼物質？
費了這麼大勁，不知講清楚了沒有。
那個文獻我也沒有，你能夠求助一下

——————————————————————————————————————————————————————

感受你講了後，確實對這個勢有了不少了解（以前特別迷糊）。我想用單層graphene和雙層graphene試着模擬Ar的nucleation過程。這個2D的如單層的a爲2.46（Angstrom），那我算的時候是否是應該就是2/（2.46*[sqrt（3）/2]*2.46）=0.38....
再好比一個fcc的（111）面，加入他的a=2，那麼計算出來的是否是應該就是4/（2^3）=0.5？這樣的話就跟選fcc的哪一個切面無關了呀（或者我仍是理解錯了）？

還有我覺得加上一個LJ9-3wall（如graphene）就是一個虛擬的wall；但若是我再在其餘位置加一個一樣組成的如graphene-wall時，而且給出這個graphene中各原子的座標，並且這個wall不是LJ9-3的勢。或者說在towhee中模擬的只要是加了一個用來吸附的表面，那麼它就應該是當作LJ9-3wall來處理嗎。這樣他的原子排列仍是很差肯定（如模擬外延生長時要選不一樣的面）？
因此加入一個本身引入的面是否是應該和加一個LJ9-3是徹底兩個概念吧？

還有，以你的經驗，你以爲像graphene有近1000做用的碳原子，這個得當作一個獨立的大分子來創建初始文件，你以爲應該用哪一種方式建呢，這樣就得在input_style下面有1000個做用的vibration手動參數輸入，感受太恐怖了。

最近在這幾個問題上實在搞不定，就只能想到麻煩你了

——————————————————————————————————————

還有一個問題是，我用ms建了graphene的模型，導出pdb文件格式，可是他是以ATOM位頭的，pdb2towhee轉換好像只能是HETATM打頭的，這個格式轉換是個問題啊。是否是應該本身編一個轉換的腳本呢。

還有就是如Ar氣在graphene上nucleation的問題，我是把Ar放在towhee_coords裏，而後把graphene放在towhee_tempalte裏面，這個石墨烯用的是trappe-ua的力場，問題是每一個碳原子都要輸入一個vibration成鍵的參數，感受工做量特別大。不知道石墨烯有沒有比較好的輸入方式。謝謝！

——————————————————————————————————————————

單層和雙層石墨烯使用9-3 Wall 勢能模擬可能不太合適，由於積分推導過程當中是認爲表面厚度無限大。1000個C原子使用全原子勢能也是能夠接受的吧。
你理解密度的計算仍是不太正確。9-3 wall 積分時認爲表面厚度無限大，那麼它的密度就等於體相的密度，跟切的表面沒有關係。
沒有明白你第二段的意思。你能再簡單地講一下嗎？
towhee中大分子的建模確實比較頭疼。不過，若是你的graphene在模擬過程當中是固定不動的，那就簡單了。你把1000個原子的石墨烯看成1000個獨立的分子，這樣就不須要輸入vibration 參數了。

————————————————————————————————————————————

pdb 轉化有問題能夠修改pdb2towhee.f90 ，也能夠編腳本實現

——————————————————————————————————

很是感謝啊，我今天下午去把這個大分子上的每一個原子當作獨立的分子，程序成功運行了。我想了下，其實每一個c原子在力場都是一個獨立體對待的，因此固定不動時把它當作不少歌分子也是可行的（下午給我樂的，呵呵）。

其實這個LJ-9-3我發現它就是一個彌散的surface吧，我那段話說的意思主要是爲了肯定這個surface表面的原子排布很差肯定，如今我想通了，它應該沒有這個原子排列的概念吧。

但願你能給我舉個例子，好比我分別用bcc與fcc的鐵的來作這個LJ-9-3 wall，加入他們的晶格常數都是2（Angtrosm），那麼這個ljfrho分別是怎麼求的，值時多少？

而ljfrho下面那個ljgsig與ljfeps的值，是否是須要從新建個模型用Fe的勢與體系的勢進行一個cycle，而後從towhee-out中看Fe與體系中各原子的sigma和epsilon值？

還有towhee mannul中關於lj-9-3介紹下面那個公式裏的sigma與epsilon應該就是ljfntypes下面的一系列ljfsig與ljfeps值吧，這樣的話，若是體系有n種與surface做用的原子，那麼在積分過程當中就應該有n中lj-9-3的積分式（參數不一樣）吧

關於那個pdb2towhee的我也搞定了，用":1, $s/ATOM  /HETATM  /g"，而後就能夠了。

————————————————————————————————————————————

有進展值得慶賀。
好比bcc的鐵，每一個晶胞的體積爲2^3=8, 包含2個Fe原子，那麼它的ljrho 應該是2/8=0.25。
ljgsig 和ljfeps 你說的是對的。
學習軟件最主要的仍是靠手冊學習，towhee的手冊寫得仍是挺詳細的，我當時也是一條條本身看的，而後結合例子進行學習。

——————————————————————————————————————————————————

關於涉及到石墨烯的建模，由於石墨烯是固定的，因此我按照你說的那種方法把每一個原子當作一個獨立的分子來建模。這樣存在着問題：
1.因爲忽然分子數多了因此對於inix， iniy， iniz 那三個值就須要增長，可是不知道這三個值究竟是什麼意思，把盒子劃分開有意義嗎（由於原子座標都是從towhee_coords中讀取的）？

2.當我用c60試驗的時候，這種方法一切正常，可是當我用涉及到1000多碳原子的石墨烯用這方法時，運行什麼的都正常，可是在輸出文件中有問題：信息只輸出到關於各原子間勢那塊就沒了，也就是顯示了LJ公式，而後是sigma   epsilon這一列，下面沒東西了。。。不知道怎麼回事，可是程序確實還在運行中，並且還不停有box_*.pdb文件產生，好生納悶啊

但願你能幫我排除下故障，謝謝啊

————————————————————————————————————————

我下午去試了下，好像是由於沒有算完，因此到那裏就不顯示了。。。

——————————————————————————————————————————————

樓主你好，我也是剛剛接觸到towhee軟件的，我是作三元萃取體系氣液平衡的，如今須要模擬純組分噻吩(thiophene)的氣液平衡狀況，在模擬以前看了一下要用的力場文件，發現裏面缺乏一部分我所須要的鍵長鍵角文件，請問那個力場文件是否是沒有顯示所有能夠處理的鍵長鍵角數據呢？（由於以前讀過一篇文獻說towhee鐘的Trappe-EH力場能夠處理N原子的問題，可是我進入到towhee 7.0.2版本的那個力場文件後發現裏面根本沒有N原子的數據。因此纔有這種疑問##）。
若是真的沒有我該本身添加力場的數據麼？須要的話是簡單地從別的聯合原子力場中找到N原子數據，而後粘貼我如今用的這個力場中嗎？
小弟初學，麻煩你了！

——————————————————————————————————————————————————

對這個力場我不是很熟悉，最好的辦法是找文獻中Trappe的參數，而後按照towhee的格式添加到力場文件中。若是Trappe力場確實沒有你所須要的參數，能夠從別的力場獲得參數，必定要二者的力場形式相同（好比都是LJ 12-6）。

————————————————————————————————————————————————————

樓主師兄你好，我又有幾個問題想問一下你。
1.兩個力場的混合規則不一樣，可不能夠用什麼方法給其中一個導入一種混合規則呢？若是能夠怎麼作？
2.對於手動輸入電荷時，若是關於力場的描述文件中的Coulombic interactions項下有電荷值，我是否是就該輸入此電荷值，若是Coulombic interactions沒有值，我是否是就不輸入電荷值呢，即爲0.0d0？入Dubb2004勢，用來創建SiO2分子時？
3,在創建一個分子時，若是有多個units的atom組成，那麼對於towhee_coords文件裏面是否是應該按照input_style下的unit值對應的把一種atom輸入完再輸入另外一種呢？
4，您知道這個程序的結果文件pdb或者movie文件裏怎麼判斷材料是fcc仍是bcc之類的嗎，或者您知道那個軟件有這個功能嗎？
5.（補充的）對於不一樣的力場文件，如argon的epsilon值都不相同，是否是意味着argon在固體狀況下，對於不一樣的力場他的晶格常數不一樣呢？
6，對於涉及同位素的，就拿H和D來講吧，若是我又H的力場文件，可不能夠把H的mass修改一下當作D的立場文件來用呢？（對於LJ勢）
很是感謝！

————————————————————————————————————————————————————————

1. 混合規則不一樣通常不推薦聯合使用，能夠測試一下看看影響是否很明顯。
2. 不考慮庫侖做用能夠把電荷寫爲0，若是所有都不考慮能夠把庫侖做用關閉
3. 不明白？？舉個例子？
4. 這個須要你本身去判斷吧？我沒有發現有這樣的功能。
5. 應該是這樣的，注意也比較一下力場的形式
6. 能夠的。

————————————————————————————————————

你好：
1.我不知道混合規則不一樣的能不能用，因此就沒用，下午能夠去試一試
2.也就是說，只要個人體系不想考了庫倫做用，我就能夠都不寫電荷值吧
3.問題解決，由於當時我導出來的sio2的pdb文件裏si和o是間隔的寫的，對於towhee須要把si的座標寫完，再寫o的座標。
4.由於我作的是氣體低溫結晶過程，因此須要判斷究竟是凍結成了什麼結構，這個很鬱悶啊
6.關於這個我知道對於vasp是能夠這麼修改的，可是對於LJ勢的力場文件，由於不知道參數和原子的質量有沒有關係，因此不知道這麼修改獲得的同位素力場可不能夠用。


問題：還想問一下，不知道您對形核過程瞭解嗎，我不知道該怎麼在模擬過程當中計算體系的形核勢壘。若是說是用成核時的體系總能量與初始的總能量相減的話，那麼這個勢壘是否是特別依賴與初始模型的原子排布（由於我是看您的教程用packmol生成的文件，因此不知道這個生成的排列對於計算勢壘時是否具備通常性）
towhee裏面有umbrella sampling和flux forwar sampling兩種抽樣嗎？若是沒有想要加的話是在編譯後的程序里加呢，仍是加了以後從新編譯呢？

很是感謝！！！！

——————————————————————————————————————————————————

6. 對於同位素：量子力學均可以使用相同的參數，分子力學固然也能夠了
對於成核不太瞭解。形核勢壘能夠經過屢次模擬把初始結構的影響下降吧。
umbrella sampling 在towhee 中是有的，後面那個不太清楚，你能夠看一下手冊。
若是要修改程序，固然須要從新編譯了！

——————————————————————————————————————

若是能夠直接把力場文件修改當作同位素來用的話，那假如我須要涉及H2和D2兩種氣體的混合，那我是否是須要把複製一個H2的力場修改後重命名或者放到其它目錄來調用呢？

————————————————————————————————————————————————————

能夠這麼作，也能夠在原來的力場中加入D原子類型

——————————————————————————————————————————————————————

你好，關於同位素力場文件，改了mass這個值就能用我有點疑問。若是說關於H和D的，質量不一樣，而後其餘什麼都相同（力場），那麼他們兩個單獨進去的模擬應該表現的狀況同樣是嗎？由於mass好像對LJ中的勢能沒有影響吧。這樣的話，其實同位素就不可分了吧？

————————————————————————————————————————————————————————

LJ 勢能確實跟質量沒有關係。不過模擬的整個過程是這個樣子的：
初始速度和座標---> 根據勢能函數求做用力 --> 根據牛頓第二定律求加速度（！！！，跟質量有關） --> 由加速度求速度 ---> 由速度求位置 --> ......循環

————————————————————————————————————————————————————————————

這個是分子動力學的吧啊，若是是蒙卡的話，那就無法區分H2與D2了嗎，是吧？

————————————————————————————————————————————————————

很差意思弄混了。

我想應該是沒有辦法區分的

——————————————————————————————————————————————

那就傻了，難道只能換lammps了嗎？我在towhee裏面的力場沒有看到有D2和T2的。。

——————————————————————————————————————————————————————

不知你研究的是同位素的動力學效應仍是熱力學效應，若是是動力學效應，MC 自己是不能作的。若是是熱力學效應，說實話，我目前還理解不了它們之間有什麼差異

——————————————————————————————————————————————————————————————————————————

我是想作H2和D2在低溫結晶過程的模擬，看他們的nucleation過程等，可是據我所知，他們的三相點是明顯不一樣的，因此氣體液化的溫度也不一樣。可是用同一個力場的話感受就是全同的了。。。體現不出什麼差異來？

——————————————————————————————————————————————————————————————————————————

你好，看來只能又麻煩你了啊。以前你給我說towhee裏面有umbrella sampling這個算法，可是我不會用它來肯定個人nucleation過程的勢壘，不知道這個具體該怎麼作，麻煩能給我說個大概的過程嗎？很是感謝，還有那個勢函數我不用去糾結了，目前先不讓作那個了，只作氬氣的。

——————————————————————————————————————————————————————————————————————————

你好，又要麻煩您了！
以前和您討論過創建石墨烯的模型，當時是把全部的C原子當作獨立的molecular來對待，即勢函數只考慮了非鍵參數，沒有考慮鍵參數。這個的前提是讓石墨烯不參與move（這樣的話把C-C做用搞成0更好，還不知道該怎麼來作，難道要單獨由幾何法獲得C-Ar參數作一個力場？）。
可是最近一個師兄說個人模型太理想了，應該讓C原子也move，若是這樣的話，就須要考慮石墨烯中的非鍵做用。若是這樣的話不知道這個石墨烯的模型該怎麼建？
若是把全部的鍵都算上的話，在towhee_input下面得輸太多太多行了。。。
關於這個問題，學長，你怎麼看？

———————————————————————————————————————————————————————————————————————————————