亭亭五月天在线观看,亭亭五月天在线观看,国产最新av一区二区,国产 高清 中文字幕,99re热久久亚洲综合精品成人,熟妇 一区二区三区,一级做a爰片性色毛片武则天,美女的骚穴视频播放,国产美女午夜免费视频

24小時(shí)熱門(mén)版塊排行榜    

查看: 16108  |  回復(fù): 63
【有獎(jiǎng)交流】積極回復(fù)本帖子,參與交流,就有機(jī)會(huì)分得作者 luyaobao 的 765 個(gè)金幣 ,回帖就立即獲得 5 個(gè)金幣,每人有 1 次機(jī)會(huì)
當(dāng)前只顯示滿足指定條件的回帖,點(diǎn)擊這里查看本話題的所有回帖

luyaobao

木蟲(chóng) (著名寫(xiě)手)


[交流] 分子動(dòng)力學(xué)模擬及其LAMMPS實(shí)現(xiàn)

(課程尚未完成,持續(xù)更新中)(三萬(wàn)字超長(zhǎng)文預(yù)警)

前言

1. 分子動(dòng)力學(xué)基礎(chǔ)

1.1 原理.

1.2 勢(shì)能

1.3 溫度、壓強(qiáng)

2. linux系統(tǒng)常見(jiàn)命令和操作

3. lammps介紹和并行版安裝(linux版)

4. lammps模擬基本流程

5. 建模和后處理軟件概述

6. 出錯(cuò)解決思路

7. 專題一:納米流動(dòng)(部分內(nèi)容)

8. 專題二:納米流體(待更新)

9. 專題三:潤(rùn)濕(部分內(nèi)容)

9.1建模

10. 專題四:拉伸、壓縮、彎曲、扭轉(zhuǎn)、剪切(待更新)

11. 專題五:傳熱(待更新)

12. 專題六:粘度和擴(kuò)散(待更新)

13. 專題八:摩擦和磨損(待更新)

14. 專題九:反應(yīng)分子動(dòng)力學(xué)reaxff(待更新)

15. 專題十:耗散粒子動(dòng)力學(xué)(待更新)

16. 專題十一:自由能(待更新)

17. 專題十二:原子沉積和鍍膜(待更新)

前言

分子動(dòng)力學(xué)模擬在科研中的使用越來(lái)越廣泛。lammps作為分子動(dòng)力學(xué)中一款優(yōu)秀的開(kāi)源軟件,使用率非常高。大量高水平的論文都是采用lammps完成的,其中不乏nature,science,prl,jacs,nature communication,pnas等頂刊。但是對(duì)于零基礎(chǔ)的同學(xué),特別是本科生和低年級(jí)碩士生,lammps的學(xué)習(xí)曲線還是比較陡峭的。很多同學(xué)都是課題組內(nèi)第一個(gè)使用分子動(dòng)力學(xué)模擬和lammps,其困難可想而知。如果能有一個(gè)針對(duì)lammps的系列教程對(duì)學(xué)生掌握這個(gè)工具會(huì)大有裨益。本教程就是希望能夠在學(xué)生使用lammps初期降低學(xué)習(xí)難度,少走彎路,快速上手,將自己的精力集中在具體的科學(xué)問(wèn)題上,而不是軟件的學(xué)習(xí)上。本教程假設(shè)你之前對(duì)linux,分子動(dòng)力學(xué)和lammps一無(wú)所知,從頭開(kāi)始講,因此比較淺顯,并未就某些問(wèn)題進(jìn)行深入介紹(其實(shí)是深入的我也不會(huì),哈哈)。

我將分子動(dòng)力學(xué)研究科學(xué)問(wèn)題分為三個(gè)部分:軟件使用,分子動(dòng)力學(xué)理論,科學(xué)問(wèn)題發(fā)現(xiàn)。如果把分子動(dòng)力學(xué)比作武林,那么軟件的使用僅僅是招式,分子動(dòng)力學(xué)理論則是內(nèi)功心法,而科學(xué)問(wèn)題則是對(duì)手。要想成為一名武林高手,在江湖上闖出一番事業(yè),既要熟練掌握各種招式,也要有很深的內(nèi)功修為,同時(shí)也要有一些值得尊敬的對(duì)手。只有用深厚的內(nèi)功催動(dòng)精妙的招式,才能打敗可敬的對(duì)手,受到江湖傳頌。本教程只是lammps的使用教程,它會(huì)教給你使用lammps的一招一式,但要將招式練得熟練,還需要自己在科研中積累大量模擬經(jīng)驗(yàn)。而內(nèi)功心法的修煉則需要長(zhǎng)時(shí)間的學(xué)習(xí)經(jīng)典教材和文獻(xiàn)。好的對(duì)手必須要廣泛閱讀文獻(xiàn),與導(dǎo)師深入溝通確定具有很好學(xué)術(shù)價(jià)值的科學(xué)問(wèn)題。那么有沒(méi)有一種方法像武俠小說(shuō)中的主角,掉下懸崖然后打通任督二脈,迅速成為高手。答案是當(dāng)然是沒(méi)有的。但是作為物理學(xué)分支,我推薦閱讀《費(fèi)曼物理講義》第一卷,可以起到類似的作用?傊,科研沒(méi)有坦途,靜心學(xué)習(xí),努力思考必有所成。

“l(fā)ammps即large-scale atomic/molecular massively parallel simulator,可以翻譯為大規(guī)模原子分子并行模擬器,主要用于分子動(dòng)力學(xué)相關(guān)的一些計(jì)算和模擬工作,一般來(lái)講,分子動(dòng)力學(xué)所涉及到的領(lǐng)域,lammps代碼也都涉及到了。lammps由美國(guó)sandia國(guó)家實(shí)驗(yàn)室開(kāi)發(fā),以gpl license發(fā)布,即開(kāi)放源代碼且可以免費(fèi)獲取使用,這意味著使用者可以根據(jù)自己的需要自行修改源代碼。lammps可以支持包括氣態(tài),液態(tài)或者固態(tài)相形態(tài)下、各種系綜下、百萬(wàn)級(jí)的原子分子體系,并提供支持多種勢(shì)函數(shù)。且lammps有良好的并行擴(kuò)展性!薄园俣劝倏。lammps的官網(wǎng)是https://www.lammps.org/。在lammps的官網(wǎng)中有大量關(guān)于軟件的信息,可以好好探索一番。有幾個(gè)地方可以經(jīng)?纯。第一是publication部分https://www.lammps.org/papers.html,記錄引用lammps的所有文章,多看看對(duì)自己的研究有所啟發(fā)。maillist部分https://www.lammps.org/mail.html是lammps的支持社區(qū)。當(dāng)你有問(wèn)題的時(shí)候可以檢索maillist,你遇到的問(wèn)題別人也遇到過(guò)。所以查看開(kāi)發(fā)者對(duì)這些問(wèn)題的解答。

如何學(xué)習(xí)lammps?你可以把本教程當(dāng)做一個(gè)開(kāi)始。lammps的手冊(cè)中詳細(xì)介紹了軟件的各個(gè)方面。一定要好好閱讀手冊(cè)。手冊(cè)的前四章要認(rèn)真閱讀。經(jīng)常使用的命令也要仔細(xì)閱讀。當(dāng)你要實(shí)現(xiàn)某種功能的時(shí)候就把手冊(cè)打開(kāi)看看命令。仔細(xì)閱讀學(xué)習(xí)lammps官方手冊(cè)是成為lammps高手的必經(jīng)之路,閱讀學(xué)習(xí)一本分子動(dòng)力學(xué)模擬經(jīng)典教材是增加修為的關(guān)鍵。希望本教程能教會(huì)你lammps的基本招式。祝你好運(yùn)!

lammps官網(wǎng)推薦了幾本分子動(dòng)力學(xué)模擬教材:

books about molecular dynamics generally or lammps specifically

note that these are not endorsements of particular books. we simply want to make the lammps user community aware of them as potentially useful resources.

general md books:

allen & tildesley - computer simulation of liquids
frenkel & smit - understanding molecular simulation: from algorithms to application
griebel, knapek, zumbusch - numerical simulation in molecular dynamics: numerics, algorithms, parallelization, applications
tuckerman - statistical mechanics: theory and molecular simulation
books about lammps specifically:

mubin & li - extending and modifying lammps (see below)

1. 分子動(dòng)力學(xué)基礎(chǔ)

1.1 原理

本節(jié)只是簡(jiǎn)要介紹分子動(dòng)力學(xué)的原理,如果要深入學(xué)習(xí)相關(guān)內(nèi)容,可閱讀前言中推薦的教材。

分子動(dòng)力學(xué)的基礎(chǔ)是牛頓力學(xué),也即經(jīng)典力學(xué)。經(jīng)典力學(xué)中有三個(gè)主要內(nèi)容:質(zhì)點(diǎn)、力和運(yùn)動(dòng)。牛頓第二定律是經(jīng)典力學(xué)的核心,下面兩個(gè)質(zhì)點(diǎn)的方程大家一定不陌生。


我們把他寫(xiě)成離散的形式,也稱為差分形式:




如果delta_t足夠小那么離散形式可以以足夠的精度近似連續(xù)形式。首先我們假設(shè)質(zhì)點(diǎn)的受力之依賴于它的位置,也就是說(shuō)只是位置的函數(shù)。事實(shí)上這是一個(gè)很普遍成立的假設(shè),比如重力或者萬(wàn)有引力,在質(zhì)量確定的情況下只依賴于兩個(gè)物體之間的相對(duì)距離。在比如靜電力,在電荷確定的情況下只依賴于兩個(gè)帶電質(zhì)點(diǎn)的相對(duì)距離。有了這個(gè)假設(shè),我們考慮一個(gè)一維的簡(jiǎn)單情形。有兩個(gè)質(zhì)點(diǎn)放置在光滑的地面上,如下圖。


在t=0時(shí)刻,質(zhì)量為m_a位于x_0^a的質(zhì)點(diǎn)a以速度v_0^a朝質(zhì)量為m_b位于x_0^b且靜止(即v_0^b=0)的質(zhì)點(diǎn)b運(yùn)動(dòng)。兩個(gè)質(zhì)點(diǎn)之間的受力f_ab只是質(zhì)點(diǎn)a與質(zhì)點(diǎn)b之間距離x_ab的函數(shù),也即只要知道x_ab,我們就能計(jì)算出來(lái)f_ab。然后,我們要問(wèn)的問(wèn)題是在后續(xù)的時(shí)間中兩個(gè)質(zhì)點(diǎn)是怎么運(yùn)動(dòng)的?我們這時(shí)設(shè)定一個(gè)很小的delta_t。我們知道初始速度那么有



這樣我們就知道了t_0+delta_t時(shí)刻,質(zhì)點(diǎn)a和b的位置和速度。由新的位置又可以得到新的受力。重復(fù)上述公式就可以知道后續(xù)所有時(shí)刻兩個(gè)質(zhì)點(diǎn)的位置和速度,也就是知道兩個(gè)質(zhì)點(diǎn)的運(yùn)動(dòng)細(xì)節(jié)?偨Y(jié)一下我們結(jié)算所有時(shí)刻質(zhì)點(diǎn)的運(yùn)動(dòng)細(xì)節(jié)要提前知道哪些信息。質(zhì)點(diǎn)的質(zhì)量、初始速度、初始位置和依賴于位置的受力函數(shù)。如果我們?cè)龠x一個(gè)足夠小的delta_t,那么我們就能夠準(zhǔn)確的得到質(zhì)點(diǎn)a和b的所有運(yùn)動(dòng)信息。

我們把情形擴(kuò)展到二維。在一個(gè)光滑的臺(tái)球桌面上,白球在初始時(shí)刻以確定的速度撞向其他擺好的臺(tái)球。如果我們把臺(tái)球都看成質(zhì)點(diǎn),又知道所有臺(tái)球的質(zhì)量、初始速度、初始位置和依賴于位置的受力函數(shù),我們就能預(yù)測(cè)擊球后所有時(shí)刻臺(tái)球的運(yùn)動(dòng)軌跡,從而判斷臺(tái)球是否能夠進(jìn)袋。我們?cè)侔亚樾螖U(kuò)展到三維?紤]宇宙中只存在太陽(yáng)和太陽(yáng)系的八大行星。我們是否能夠通過(guò)萬(wàn)有引力去預(yù)測(cè)所有行星的軌跡。答案當(dāng)然是可以的。只是此時(shí)我們知道某個(gè)行星除了來(lái)自太陽(yáng)的引力,還有來(lái)自其他行星的應(yīng)力。我們?cè)谟?jì)算某個(gè)行星的受力時(shí),要考慮該行星與其余所有天體之間的受力。這就是經(jīng)典力學(xué)能夠告訴我們的東西。

在這個(gè)概念上,拉普拉斯提出了一個(gè)想法。如果存在一只神奇的動(dòng)物,它能基于經(jīng)典力學(xué)解算出世界上所有東西的位置和速度并能加以分析,那么世界上發(fā)生的一切及其背后的原因都一覽無(wú)余的呈現(xiàn)在它眼前。這就是經(jīng)典力學(xué)的世界觀和認(rèn)識(shí)論。分子動(dòng)力學(xué)的理論就是基于這樣的認(rèn)識(shí)方法。像上面的例子一樣,分子動(dòng)力學(xué)的目的就是為了求解任意時(shí)刻所有組成物體原子的位置和速度,從而認(rèn)識(shí)物質(zhì)世界。


為什么這個(gè)方法具有很大的威力呢?這是多年科學(xué)實(shí)踐告訴我們的。事實(shí)上著名物理學(xué)家理查德·費(fèi)曼對(duì)這個(gè)問(wèn)題有著精彩的描述。理查德·費(fèi)曼在其傳世物理教材《費(fèi)曼物理講義》中提出這樣一個(gè)有意思的問(wèn)題:“假如由于某種大災(zāi)難,所有的科學(xué)知識(shí)都丟失了,只有一句話傳給下一代,那么怎樣才能用最少的詞匯來(lái)表達(dá)最多的信息呢?”費(fèi)曼認(rèn)為這句話應(yīng)該是原子假設(shè):“所有物體都是用原子構(gòu)成的——這些原子是一些小小的粒子,它們一直不停地運(yùn)動(dòng)著,當(dāng)彼此略微離開(kāi)地時(shí)互相吸引,當(dāng)彼此過(guò)于靠近時(shí)又相互排斥!边@句話之所以重要是因?yàn)榛谶@個(gè)原子假設(shè)我們能解釋巨量自然界中的現(xiàn)象。比如,氣體壓強(qiáng)來(lái)自于大量氣體分子對(duì)壁面撞擊;水蒸發(fā)吸熱來(lái)自于水分子從體相水中逃逸至外界環(huán)境而帶走能量所致。。。。。更多詳細(xì)例子可閱讀《費(fèi)曼物理講義》第一卷的相關(guān)內(nèi)容。

費(fèi)曼把原子假設(shè)放在如此高度的原因總結(jié)起來(lái)就是——掌握組成物體內(nèi)部原子的行為就可解釋和理解物體的性質(zhì)和行為。只是這個(gè)時(shí)候我們關(guān)注的是大量原子在一起表現(xiàn)出來(lái)的整體特性,而不是關(guān)注具體某個(gè)原子的運(yùn)動(dòng)細(xì)節(jié)。事實(shí)上,實(shí)踐告訴我們大量原子在一起就會(huì)表現(xiàn)出某些特定的確定的性質(zhì)和行為,也就我們實(shí)際生活生產(chǎn)中能直接測(cè)量和感受的性質(zhì)和行為。到這里我們可以看到分子動(dòng)力學(xué)實(shí)際上就是一只拉普拉斯獸。分子動(dòng)力學(xué)的基本任務(wù)就是獲取物體在任意時(shí)刻組成原子的所有位置和動(dòng)量然后利用統(tǒng)計(jì)力學(xué)知識(shí)理解物體的性質(zhì)和行為。如果這樣還不能說(shuō)服你分子動(dòng)力學(xué)的實(shí)力。就講一個(gè)事實(shí)——2013年諾貝爾化學(xué)獎(jiǎng)?lì)C給了創(chuàng)立跨尺度分子模擬的三位科學(xué)家。通過(guò)分子動(dòng)力學(xué)模擬我們可以深刻理解我們所研究的對(duì)象。目前各種硬件超強(qiáng)的計(jì)算能力為分子動(dòng)力學(xué)創(chuàng)建了一個(gè)大有可為的時(shí)代。

分子動(dòng)力學(xué)的基本任務(wù)就是獲得研究對(duì)象不同時(shí)刻的位置和動(dòng)量,然后基于統(tǒng)計(jì)力學(xué)知識(shí)獲得想要的物理量,解釋對(duì)象的性質(zhì)和行為。因此,分子動(dòng)力學(xué)的模擬流程超級(jí)簡(jiǎn)單。第一步,設(shè)置研究對(duì)象組成原子或者粒子(下面統(tǒng)一稱為粒子)初始位置和速度;第二步,基于粒子的位置計(jì)算每個(gè)粒子受到的合力,并基于牛二定律計(jì)算粒子的加速度;第三步,基于加速度計(jì)算粒子下一時(shí)刻的速度;第四步,基于下一時(shí)刻的速度計(jì)算下一時(shí)刻的位置;第五步,循環(huán)第二步到第四步的過(guò)程,得到一系列時(shí)刻粒子的位置和速度;第六步,基于位置和速度信息得到描述對(duì)象性質(zhì)和行為的物理量。


分子動(dòng)力學(xué)的基本原理就是這么簡(jiǎn)單,正因?yàn)楹?jiǎn)單所以有效。需要指出的是在實(shí)際編程求解的過(guò)程考慮到計(jì)算精度等原因上述過(guò)程在某些細(xì)節(jié)上會(huì)略有改動(dòng)。但是,基本上你可以認(rèn)為這就是分子動(dòng)力學(xué)的原理和流程。上述過(guò)程由lammps幫你實(shí)現(xiàn),你需要做的就是設(shè)置研究對(duì)象的基本信息,告訴給lammps,然后開(kāi)始模擬。

1.2 勢(shì)能

前面我們知道要獲得所有原子的軌跡,需要知道原子的質(zhì)量,初始位置,初始速度和依賴于位置的受力函數(shù)。原子的質(zhì)量,初始速度和初始位置都是我們?nèi)藶樵O(shè)定的值。不同研究對(duì)象的原子質(zhì)量,初始速度和位置很容易獲得。不同研究對(duì)象更本質(zhì)的差別在于依賴于位置的受力函數(shù),這個(gè)函數(shù)又被稱為相互作用。相互作用的形式和參數(shù)決定了研究對(duì)象的性質(zhì)和行為。在分子動(dòng)力學(xué)中,這個(gè)相互作用通過(guò)原子之間的勢(shì)能進(jìn)行描述,因此又稱為勢(shì)能函數(shù)或勢(shì)函數(shù)。而受力是勢(shì)能對(duì)位置的負(fù)梯度,因?yàn)閷?shí)踐表明原子之間的受力都是保守力。舉個(gè)例子,重力和重力勢(shì)能。取豎直向上為正方向,則重力表達(dá)式為


重力勢(shì)能的表達(dá)式是


即重力勢(shì)能對(duì)位置的負(fù)梯度就是重力


為什么要用勢(shì)能,而不是受力來(lái)表達(dá)相互作用呢?這是因?yàn)樵诮^大多數(shù)模擬中采用的是笛卡爾坐標(biāo)。勢(shì)能是標(biāo)量,而力是矢量,數(shù)學(xué)上處理起來(lái)更容易。在計(jì)算受力的時(shí)候,直接用勢(shì)能對(duì)兩個(gè)原子之間的相對(duì)位置進(jìn)行求導(dǎo),然后將這個(gè)力乘以方向余弦就能得到三個(gè)方向上的受力。舉例說(shuō)明,在三維空中有兩個(gè)原子


兩個(gè)原子之間的距離


那么原子對(duì)原子的受力就為


這樣不管是程序的簡(jiǎn)潔性還是計(jì)算效率都會(huì)很好。所以要開(kāi)展分子動(dòng)力學(xué)計(jì)算,首先必須確定所有原子之間的勢(shì)能函數(shù),即勢(shì)能關(guān)于兩個(gè)原子之間相對(duì)位置的函數(shù)。確定函數(shù)分為兩個(gè)部分:函數(shù)表達(dá)式和表達(dá)式中的常數(shù)值。對(duì)于一個(gè)對(duì)象,描述其相互作用的勢(shì)能函數(shù)表達(dá)式和表達(dá)式中常數(shù)值的集合稱為該對(duì)象的力場(chǎng)。通常一種力場(chǎng)可以用來(lái)描述一類對(duì)象。分子動(dòng)力學(xué)模擬一般研究的對(duì)象可以分為:惰性氣體,分子物質(zhì),金屬晶體和非金屬晶體。分子動(dòng)力學(xué)模擬絕大多數(shù)情況用來(lái)研究固體和液體,即凝聚態(tài)。對(duì)于氣體分子動(dòng)力學(xué)模擬效率太低,有更好的辦法進(jìn)行研究。簡(jiǎn)單液體原子之間相互作用只具有范德瓦爾斯作用的物質(zhì)。描述范德瓦爾斯作用最常用的表達(dá)式是lennard-jones12/6勢(shì)能。通常用lj12/6描述的物質(zhì)也稱為簡(jiǎn)單液體。





上圖是lj/126勢(shì)函數(shù)的圖像。對(duì)于該勢(shì)函數(shù)有幾個(gè)重要的點(diǎn)要說(shuō)一下。勢(shì)能的最低點(diǎn)的到零點(diǎn)之間的絕對(duì)值為ε,稱為勢(shì)能函數(shù)的特征能量,勢(shì)能的零點(diǎn)位于處r=σ,稱為勢(shì)能函數(shù)的特征長(zhǎng)度。勢(shì)能函數(shù)的最低點(diǎn)位于處r=2^1/6σ,在該點(diǎn)處勢(shì)能對(duì)位置的導(dǎo)數(shù)是零,所以該點(diǎn)處受力也為零。該點(diǎn)也稱為平衡位置,我們可以看到在平衡位置出勢(shì)能是負(fù)的。所以當(dāng)體系處于平衡狀態(tài)時(shí),原子之間的相對(duì)距離大概也在平衡位置附近,所以勢(shì)能為負(fù)的。這也就是為什么大多數(shù)情況分子動(dòng)力學(xué)算出來(lái)的勢(shì)能都是負(fù)的。在r<2^1/6σ,兩個(gè)原子之間表現(xiàn)為相互排斥,在r>2^1/6σ,兩個(gè)原子之間表現(xiàn)為相互吸引。從勢(shì)函數(shù)的圖像上可以看出在較大的時(shí)候,勢(shì)能趨于零,其導(dǎo)數(shù)也趨于零。也就是當(dāng)兩個(gè)原子距離很遠(yuǎn)的時(shí)候,原子之間的相互作用可以忽略不計(jì)。因此在實(shí)際計(jì)算的時(shí)候都會(huì)使用截?cái)嗟膌j/126。截?cái)嗑嚯x成為截?cái)喟霃剑渲荡蟾旁。就像重力?shì)能一樣,勢(shì)能的計(jì)算需要一個(gè)參考位置。lj/126的參考位置選擇在了無(wú)限遠(yuǎn)處。類似lj/126勢(shì)能以及其他拓展的lj勢(shì)能,都稱為對(duì)勢(shì),其含義是兩個(gè)原子之間的相互作用只由這一對(duì)原子決定,與其它相鄰原子沒(méi)有關(guān)系。分子物質(zhì)原子之間有化學(xué)鍵連接的物質(zhì)。比如常見(jiàn)的水,甲烷,乙醇,醋酸,氯化鈉,聚乙烯,酯類,氨基酸,rna,dna等。分子物質(zhì)是分子動(dòng)力學(xué)模擬應(yīng)用最廣的對(duì)象之一。為了描述分子物質(zhì)研究者開(kāi)發(fā)出了一系列力場(chǎng)。值得說(shuō)明的是力場(chǎng)是一個(gè)經(jīng)驗(yàn)性的東西。研究者通過(guò)自己的物理直覺(jué)和反復(fù)試驗(yàn),最終建立起適合某一類物質(zhì)的力場(chǎng)。對(duì)于分子物質(zhì),相互作用包括兩類:非鍵結(jié)勢(shì)能和鍵結(jié)勢(shì)能。非鍵結(jié)用來(lái)描述沒(méi)有化學(xué)鍵連接原子之間的相互作用,鍵結(jié)勢(shì)能用來(lái)描述有化學(xué)鍵連接原子之間的相互作用,如下圖。


在一個(gè)體系中任意兩個(gè)非鍵結(jié)的原子,如果其距離小于截?cái)喟霃侥敲炊家紤]它們之間的非鍵結(jié)相互作用。對(duì)于大多數(shù)力場(chǎng)非鍵結(jié)勢(shì)能用lj12/6描述。對(duì)于鍵結(jié)原子之間的相互作用則通過(guò)鍵結(jié)勢(shì)函數(shù)進(jìn)行描述。鍵結(jié)相互作用分為鍵能(bond),鍵角能(angle),二面角能(dihedral),離平面能(improper)。兩個(gè)以化學(xué)鍵相連的兩個(gè)原子都具有鍵能相互作用,當(dāng)化學(xué)鍵偏離平衡長(zhǎng)度的時(shí)候兩個(gè)原子就會(huì)出現(xiàn)相互作用。三個(gè)化學(xué)鍵連接的相鄰原子之間均具有鍵角能,當(dāng)三個(gè)原子之間的角度偏離平衡角度的時(shí)候就會(huì)產(chǎn)生鍵角能。四個(gè)化學(xué)鍵連接的相鄰原子之間會(huì)存在二面角能相互作用,當(dāng)四個(gè)原子構(gòu)成的二面角偏離平衡二面角的時(shí)候四個(gè)原子之間會(huì)出現(xiàn)二面角相互作用。四個(gè)相鄰化學(xué)鍵連接的相鄰原子之間有時(shí)還會(huì)存在離平面勢(shì)能作用,當(dāng)?shù)谒膫(gè)原子偏離前三個(gè)原子所構(gòu)成平面的平均夾角時(shí),四個(gè)原子之間會(huì)存在離平面能相互作用。不同的分子會(huì)具有不同的鍵結(jié)勢(shì)能。如氧氣分子就只具有鍵能相互作用,水分子具有鍵能和鍵角相互作用,十二烷分子具有鍵能,鍵角能和二面角相互作用分子。元素構(gòu)成較為復(fù)雜的分子,如rna,就會(huì)同時(shí)具有四種鍵結(jié)相互作用。具體一個(gè)分子具有哪些鍵結(jié)能相互作用,由描述該分子的力場(chǎng)所定義。你只需要使用合適的軟件在建模的時(shí)候會(huì)自動(dòng)定義好各種鍵結(jié)信息和勢(shì)能函數(shù)的參數(shù)。值得指出的是,在分子體系中原子一般都帶有電荷,所以非鍵結(jié)原子之間會(huì)存在靜電相互作用。在大多數(shù)力場(chǎng)中,原子的電荷都是部分電荷,也就是電荷一般都是零點(diǎn)幾幾幾。這是由于在共價(jià)鍵中任何一個(gè)原子都并沒(méi)有完全得到或失去一個(gè)電子,只是部分的得到或失去電子,因此在描述其電荷的時(shí)候都是零點(diǎn)幾幾幾。少數(shù)強(qiáng)電解質(zhì)會(huì)有整數(shù)電荷,比如氯化鈉中的氯離子和鈉離子,分別是-1和+1電荷。這里的電荷單位是一個(gè)電子或質(zhì)子所帶的電荷量,也就是單位電荷量。

為什么會(huì)采用多種鍵結(jié)勢(shì)能函數(shù)來(lái)描述分子體系呢?這是因?yàn)殒I結(jié)的原子之間電子云的相互作用會(huì)對(duì)分子的構(gòu)型進(jìn)行限制,導(dǎo)致分子會(huì)具有特定的幾何形狀。分子動(dòng)力學(xué)為描述這種情形就采用這樣鍵能組合的形式。這是一種經(jīng)驗(yàn)的方法。實(shí)踐表明這是一種很好的辦法。只要有合適的勢(shì)函數(shù)形式和參數(shù),采用的力場(chǎng)就能夠準(zhǔn)確描述分子的構(gòu)型和之間的相互作用。目前研究者針對(duì)分子體系已經(jīng)開(kāi)發(fā)出多套力場(chǎng)。
金屬晶體和非金屬晶體

這兩類物質(zhì)具有特殊的結(jié)構(gòu),研究者為它們開(kāi)發(fā)了額外的勢(shì)函數(shù),稱為多體勢(shì)函數(shù),來(lái)描述原子之間的作用。所謂多體就是說(shuō)兩個(gè)原子之間的相互作用不僅取決于這兩個(gè)原子還與其周圍的原子有關(guān)系。針對(duì)某種物質(zhì)研究者已經(jīng)開(kāi)發(fā)出相應(yīng)的力場(chǎng)參數(shù),只需要搜索一般就能找到對(duì)應(yīng)的勢(shì)函數(shù)參數(shù)或文件。我對(duì)這些勢(shì)函數(shù)不太了解,就不在這里贅述了。

以上就是關(guān)于勢(shì)能和勢(shì)函數(shù)的內(nèi)容。需要指出的是為了提高模擬精度,研究者還針對(duì)開(kāi)發(fā)出了針對(duì)某一種物質(zhì)的勢(shì)函數(shù)。當(dāng)你要模擬這些物質(zhì)時(shí),直接查找原始文件中的參數(shù)就行。關(guān)于力場(chǎng)的選擇,當(dāng)然是選針對(duì)要研究體系的力場(chǎng)。但是對(duì)于lammps來(lái)說(shuō)建模是個(gè)很麻煩的事。如果你能找到一種力場(chǎng)描述的體系,并且成功建模,那你就偷著樂(lè)吧。不同力場(chǎng)之間的差別并不是那么明顯。事實(shí)上,分子動(dòng)力學(xué)模擬最大的價(jià)值就在于給出研究對(duì)象的一般性規(guī)律和背后的機(jī)理,而不是計(jì)算出某些精確的數(shù)值。而大多數(shù)力場(chǎng)能夠捕捉研究對(duì)象的本質(zhì)特征,這就基本上能夠保證這個(gè)力場(chǎng)能夠再現(xiàn)研究對(duì)象的本質(zhì)規(guī)律,這對(duì)于我們來(lái)說(shuō)就已經(jīng)足夠了。當(dāng)然,力場(chǎng)的開(kāi)發(fā)一直在發(fā)展中以越來(lái)越精確的描述研究對(duì)象。比如幾年來(lái)發(fā)展出來(lái)的機(jī)器學(xué)習(xí)力場(chǎng)就能夠以量子的精度開(kāi)展分子動(dòng)力學(xué)模擬。這是分子動(dòng)力學(xué)模擬方法中的一大進(jìn)展。當(dāng)你要驗(yàn)證一個(gè)力場(chǎng)的時(shí)候,計(jì)算研究對(duì)象你所關(guān)心的物理量的數(shù)值和實(shí)驗(yàn)或者文獻(xiàn)對(duì)比,如果接近,力場(chǎng)就能用。

1.3 溫度、壓強(qiáng)(這里寫(xiě)的有點(diǎn)簡(jiǎn)單了,先就這樣吧) 溫度對(duì)于分子動(dòng)力學(xué)模擬是一個(gè)很重要的物理量,因?yàn)榉肿觿?dòng)力學(xué)一般用來(lái)模擬凝聚態(tài),如果溫度過(guò)高就氣化了。初中物理告訴我們溫度是原子熱運(yùn)動(dòng)劇烈程度的度量。原子熱運(yùn)動(dòng)越劇烈,溫度越高。也就是說(shuō)溫度和原子的熱運(yùn)動(dòng)速度是直接相關(guān)的。根據(jù)統(tǒng)計(jì)力學(xué)三維系統(tǒng)中溫度和原子速度具有以下對(duì)應(yīng)關(guān)系


其中,是n原子總個(gè)數(shù),k_b是玻爾茲曼常數(shù)。通過(guò)調(diào)整原子的速度就可以實(shí)現(xiàn)控制系統(tǒng)的溫度。

壓強(qiáng)是另一個(gè)重要的物理量。壓強(qiáng)在分子動(dòng)力學(xué)中這樣定義的


其中,v是系統(tǒng)體積。第一項(xiàng)來(lái)自原子熱運(yùn)動(dòng)的貢獻(xiàn)。第二項(xiàng)來(lái)自原子之間相互作用的貢獻(xiàn),稱為維里項(xiàng)。通過(guò)調(diào)整系統(tǒng)的體積和原子的位置就可以實(shí)現(xiàn)控制系統(tǒng)的壓強(qiáng)。對(duì)于液體來(lái)說(shuō)壓強(qiáng)就等于我們平時(shí)所接觸到的壓強(qiáng)。而對(duì)于固體來(lái)說(shuō)這里的壓強(qiáng)實(shí)際上就是內(nèi)應(yīng)力。事實(shí)上,壓強(qiáng)其實(shí)就是系統(tǒng)內(nèi)部應(yīng)力張量的對(duì)角線分量。壓強(qiáng)和溫度都是強(qiáng)度量,就是說(shuō)其值與系統(tǒng)規(guī)模無(wú)關(guān)。與之相對(duì)應(yīng)的是廣延量,比如能量,當(dāng)系統(tǒng)規(guī)模增大時(shí),系統(tǒng)能量增加。

更多分子動(dòng)力學(xué)的細(xì)節(jié)將結(jié)合lammps使用過(guò)程進(jìn)行介紹。

2. linux系統(tǒng)常見(jiàn)命令和操作

lammps是在針對(duì)linux系統(tǒng)開(kāi)發(fā)的。雖然windows系統(tǒng)也能用但是會(huì)出現(xiàn)很多問(wèn)題,所以強(qiáng)烈推薦你使用linux系統(tǒng)運(yùn)行l(wèi)ammps。要搞模擬你一定要熟悉linux的命令行操作邏輯,要習(xí)慣使用腳本來(lái)寫(xiě)lammps命令。很多軟件包括開(kāi)源的和商業(yè)的都是采用腳本控制的方式來(lái)執(zhí)行的。

在linux系統(tǒng)中使用最多的不是鼠標(biāo)點(diǎn)點(diǎn)點(diǎn),而是在終端中輸入命令行來(lái)執(zhí)行操作。執(zhí)行一條命令是指在終端中輸入該命令,然后回車執(zhí)行。在ubuntu系統(tǒng)中可以使用ctrl+alt+t來(lái)打開(kāi)終端,其他linux操作系統(tǒng)可以在軟件列表中點(diǎn)擊打開(kāi)終端。這里只列舉最常用的命令。

cd命令用來(lái)進(jìn)入某個(gè)目錄,比如cd /home就是進(jìn)入home文件夾。在linux中.表示當(dāng)前目錄,..表示上一級(jí)目錄。ls命令列出當(dāng)前目錄下的所有文件和文件夾。pwd命令返回當(dāng)前目錄的路徑。mkdir創(chuàng)建文件夾,如mkdir software表示在當(dāng)前目錄下創(chuàng)建software文件夾,mkdir /home/softeare,表示在/home文件夾下創(chuàng)建software文件夾。rm命令刪除文件或文件夾。tar解壓縮命令當(dāng)壓縮文件的擴(kuò)展名是.tar.gz時(shí)可以使用tar -xvzf 壓縮文件名,將壓縮文件解壓。which命令,查找可執(zhí)行文件的路徑,如which lmp_g++_mpich就是查找可執(zhí)行文件lmp_g++_mpich的所在路徑,如果返回路徑那么這個(gè)命令就可以被系統(tǒng)找到,如果什么也沒(méi)返回就是找不到。echo命令返回變量的值,如echo $path就是返回path的值。環(huán)境變量是指系統(tǒng)或軟件定義的變量,一般用于文件查找。最常見(jiàn)的環(huán)境變量就是path,該變量記錄了眾多可執(zhí)行文件的位置。如果一個(gè)可執(zhí)行文件的路徑被添加到path中了,那么在運(yùn)行該文件時(shí)可以直接在終端中寫(xiě)該文件的名字,如果沒(méi)有被添加那么就要寫(xiě)出該文件的絕對(duì)路徑才能執(zhí)行該文件。將一個(gè)路徑添加到path中需要修改.bashrc文件。輸入執(zhí)行vim ~/.bashrc,就可以打開(kāi).bashrc文件,按i鍵進(jìn)入編輯模式,將光標(biāo)定位在文件末尾并添加以下一行

path=$path:<路徑名>,按esc鍵退出編輯模式,輸入:wq后回車,這樣就保存好了剛才的編輯,并退出了.bashrc文件回到了終端。在終端輸入source ~/.bashrc就可以是剛才的操作發(fā)生作用。sudo是臨時(shí)取得root權(quán)限命令,當(dāng)你要在系統(tǒng)目錄中添加或修改文件的時(shí)候就需要把sudo加到執(zhí)行命令前面,如sudo cp ./lmp /usr/bin就是把當(dāng)前目錄下的lmp文件拷貝到系統(tǒng)目錄/usr/bin中。帶有sudo的命令回車執(zhí)行后需要輸入用戶密碼,當(dāng)你在輸入密碼的時(shí)候終端上不會(huì)顯示密碼,但是確實(shí)是輸進(jìn)去了,輸完回車如果密碼正確就開(kāi)始執(zhí)行命令了。其他命令可以百度能夠找到大量資料。當(dāng)你在使用linux系統(tǒng)時(shí)需要問(wèn)一個(gè)問(wèn)題,你就百度,一定能找到答案。

3. lammps介紹和并行版安裝(linux版)

這里介紹在ubuntu在編譯并行版lammps并執(zhí)行。首先假設(shè)你已經(jīng)有了一臺(tái)安裝了ubuntu系統(tǒng)的電腦。你的用戶名叫me。電腦聯(lián)網(wǎng)。分別在lammps官網(wǎng),mpich2官網(wǎng)和fftw官網(wǎng)下載對(duì)應(yīng)的lammps,mpich和fftw的最新版源代碼,然后

第一步按ctrl+alt+t打開(kāi)一個(gè)終端,輸入執(zhí)行以下命令

sudo apt-get update

sudo apt-get install build-essential

sudo apt-get install vim

第二步輸入執(zhí)行

mkdir /home/me/software

mkdir /home/me/software/mpich

mkdir /home/me/software/fftw

第三步將下載好的三個(gè)源代碼的壓縮文件拷貝進(jìn)/home/me/software目錄下

第四步輸入執(zhí)行

tar -xvf <lammps壓縮文件的名稱>

tar -xvf <mpich壓縮文件的名稱>

tar -xvf <fftw壓縮文件的名稱>

第五步進(jìn)入解壓出來(lái)的mpich的目錄執(zhí)行

./configure -prefix=/home/me/software/mpich --disable-fortran

make

make install

第六步進(jìn)入解壓出來(lái)的fftw的目錄執(zhí)行

./configure -prefix=/home/me/software/fftw

make

make install

第七步輸入執(zhí)行

vim ~/.bashrc

打開(kāi).bashrc文件后按i鍵進(jìn)入邊界模式,將光標(biāo)定位在文件末尾并添加以下一行

path=$path:/home/me/software/mpich/bin

按esc鍵退出編輯模式,輸入:wq后回車,這樣就保存好了剛才的編輯,并退出了.bashrc文件回到了終端。在終端輸入

source ~/.bashrc

which mpirun

若輸出了mpirun的目錄則mpich安裝配置成功

第八步終端輸入執(zhí)行以下命令

cd /home/me/software/解壓出來(lái)lammps文件夾名稱/src

第九步打開(kāi)目錄/home/me/software/解壓出來(lái)lammps文件夾名稱/src/make/options下的makefile.g++_mpich文件。在文件以下位置修改為

mpi_inc = -dmpich_skip_mpicxx -dompi_skip_mpicxx=1 -i/home/me/software/mpich/include

mpi_path = -l/home/me/software/mpich/lib

mpi_lib =

fft_inc = -i/home/me/software/fftw/include

fft_path = -l/home/me/software/fftw/lib

fft_lib =

保存退出。

第十步在終端依次執(zhí)行

make yes-molecule

make yes-kspace

make yes-rigid

make yes-manybody

make g++_mpich

等待編譯完成src目錄下會(huì)生成lmp_g++_mpich的可執(zhí)行文件,將此可執(zhí)行文件拷貝至你的in文件所在的位置輸入

mpirun -np <并行核數(shù)> ./lmp_g++_mpich -in <in文件名字>

運(yùn)行正常則編譯成功。

4. lammps模擬基本流程

lammps的模擬需要寫(xiě)一個(gè)腳本,稱之為in文件。我將一個(gè)典型的in文件書(shū)寫(xiě)分為六個(gè)部分:1、定義模擬的基本規(guī)則,2、定義對(duì)象的幾何模型,3、定義對(duì)象的力場(chǎng)參數(shù),4、定義分組信息,5、弛豫部分,6、正式模擬過(guò)程。

第一部分 定義模擬的基本規(guī)則

根據(jù)前面說(shuō)的分子動(dòng)力學(xué)原理,那么第一步就是設(shè)置粒子的初始位置和速度?等等。在這之前我們還需要做一些前期準(zhǔn)備工作。首先是設(shè)置單位制,即涉及到的物理量的單位都是什么。為什么要這么做呢?就是為了能夠定量描述我們的對(duì)象;叵朐诟咧械闹麑(shí)驗(yàn)——小球平拋。我們?cè)趺疵枋鲂∏虻倪\(yùn)動(dòng)呢?我們會(huì)說(shuō)小球經(jīng)過(guò)多少秒,向前運(yùn)動(dòng)了多少米,下落了多少米。這里就采用了國(guó)際制單位描述小球的運(yùn)動(dòng)。只有先確定了單位制,我們才好定量描述研究對(duì)象。不同的研究對(duì)象,為了描述方便會(huì)采用不同的單位制。比如描述高鐵的速度我們會(huì)用km/h,而不是m/s。同樣,在lammps中面對(duì)不同的研究對(duì)象,lammps內(nèi)置了不同的單位制。所以in文件中的

第一條命令要寫(xiě)的就是設(shè)置單位制,即寫(xiě)units命令。語(yǔ)法就是、

units <style> #舉例:units real

lammps所有的單位制可在手冊(cè)中查到。常見(jiàn)的單位系統(tǒng)有,描述分子體系的real,描述金屬體系的metal,描述簡(jiǎn)單液體的lj。

第二條命令是設(shè)置邊模擬維度,即你模擬的是三維問(wèn)題還是二維問(wèn)題。通常我們的模擬都是三維的。語(yǔ)法是

dimension <n> #舉例:dimension 3

第三條命令是設(shè)置邊界條件,即boundary命令。通常我們希望描述對(duì)象的體相性質(zhì)。但是,看看阿佛加德羅常數(shù)的值我們就可以知道實(shí)際的對(duì)象都會(huì)包含超級(jí)多的粒子。這么多的粒子的愛(ài)計(jì)算機(jī)無(wú)法承受。那怎么辦呢?分子動(dòng)力學(xué)中采用了有限代替無(wú)限的做法。就是取對(duì)象的一部分然后加上周期性邊界條件來(lái)代替幾乎無(wú)限大的對(duì)象。如圖中所示所謂周期性邊界條件就是粒子如果從邊界穿出去,就會(huì)從對(duì)面穿回來(lái)。這樣就實(shí)現(xiàn)了以有限的粒子數(shù)代替了無(wú)限大的對(duì)象。由邊界所圍起來(lái)的區(qū)域通常稱為模擬盒子。模擬盒子的形狀最常用的就是長(zhǎng)方體。之所以采用長(zhǎng)方體是因?yàn)榫幊倘菀滋幚怼?br />

設(shè)置邊界條件在分子動(dòng)力學(xué)模擬中優(yōu)先級(jí)很高,所以第二條命令要寫(xiě)邊界條件,語(yǔ)法是

boundary <x> <y> <z> #舉例:boundary p p p

boundary p p p就是指三個(gè)方向上都是周期性邊界條件。這個(gè)設(shè)置在lammps中最常用。lammps中還提供了其他類型的邊界條件,常用的有固定邊界條件(f)和縮放邊界條件(s)。當(dāng)我們需要在某個(gè)方向上不想采用周期性時(shí)就可用f或s。f作用就是粒子從邊界跑出去就回不來(lái)了,計(jì)算就把它丟掉了。s的作用是邊界的位置會(huì)隨著粒子的位置而縮放,當(dāng)最靠近邊界的粒子往外移動(dòng)時(shí)邊界位置也跟著移動(dòng),反之亦然。需要指出的是每個(gè)方向上都有兩個(gè)邊界,當(dāng)寫(xiě)一個(gè)參數(shù)是如p就表示該方向上兩個(gè)邊界都是周期性,f和s也如此。但是,還可以這么寫(xiě)boundary p fs p,這就表示在x和z方向上是周期性邊界,在y方向上下面是f,上面是s。因此寫(xiě)該條命令時(shí)根據(jù)需要進(jìn)行設(shè)置。

第四條命令是設(shè)置建立鄰居列表的,即neighbor命令。根據(jù)分子動(dòng)力學(xué)原理我們需要計(jì)算每個(gè)原子的受力。而任何一個(gè)粒子的受力都來(lái)自于其周圍粒子給它的合力。但是,我們知道只有距離較近的粒子才會(huì)對(duì)該粒子產(chǎn)生受力,距離較遠(yuǎn)的粒子給與的力可以忽略不計(jì)。這就是鄰居的含義。所以要計(jì)算一個(gè)粒子所受的合力,首先要找到該粒子周圍有哪些鄰居,即建立鄰居列表。每一個(gè)粒子都有自己的鄰居列表。neighbor的語(yǔ)法是

neighbor <skin> <style> #舉例:neighbor 3 bin

所謂鄰居就是距離中心粒子一定范圍內(nèi)的粒子,這個(gè)范圍稱之為作用力的截?cái)喟霃剑╮_c)。也就是以中心粒子為球心,r_c為半徑內(nèi)的所有粒子都是鄰居。skin是比r_c多出來(lái)的一部分,所以建立鄰居列表的時(shí)候是以r_c+skin為半徑進(jìn)行尋找的。為什么要這么做呢?這是因?yàn)榻⑧従恿斜硎峭M(fèi)時(shí)間的一件事,那么好不容易建立起來(lái)那就盡可能多用幾次。而事實(shí)上由于分子動(dòng)力學(xué)模擬的時(shí)間步長(zhǎng)較短,粒子在一個(gè)時(shí)間步長(zhǎng)內(nèi)走不了多少距離。那么只要鄰居列表中的原子移動(dòng)距離都小于skin,那么這個(gè)列表就可以在下一步中繼續(xù)使用。這就是skin的含義。后面style指的是建立鄰居列表的方法,最常用的就是bin方法。這個(gè)方法的具體含義可以查看分子動(dòng)力學(xué)教材。入門(mén)只需要知道bin在大多數(shù)情況下都是最快的,用它就沒(méi)錯(cuò)了。而r_c在后面的命令中設(shè)置。


第五條命令是設(shè)置鄰居列表的更新頻率,即neigh_modify命令。lammps提供了很多方式設(shè)置鄰居列表的更新頻率。如果鄰居中的粒子已經(jīng)不在r_c+skin范圍內(nèi)了,那么鄰居列表就必須要更新了。當(dāng)需要更新而未更新時(shí)就會(huì)出現(xiàn)dangerous builds。當(dāng)算完之后lammps輸出的dangerous builds不為零,那么你就要修改neigh_modify命令了,或者加大skin的值。這條命令使用比較復(fù)雜,需要根據(jù)模擬經(jīng)驗(yàn)和具體的模擬測(cè)試確定最優(yōu)值。對(duì)于新手,可以直接使用以下設(shè)置:

neigh_modify every 1 delay 0 check yes

以上都是計(jì)算設(shè)置,每個(gè)模擬這些設(shè)置都大同小異。下面開(kāi)始就是針對(duì)模擬對(duì)象的設(shè)置。

第六條命令設(shè)置原子類型,即atom_style命令。我們模擬的基本構(gòu)成是原子。不同的對(duì)象粒子所具有的屬性是不同的。對(duì)于某個(gè)對(duì)象,粒子不需要的屬性就不需要設(shè)置,以節(jié)省內(nèi)存加速計(jì)算。比如研究對(duì)象是分子(如水分子),那么構(gòu)成分子的粒子(style為full)就需要具有以下屬性:粒子編號(hào),粒子所屬分子的編號(hào),粒子類別(atom type),帶電量,x坐標(biāo),y坐標(biāo),z坐標(biāo)。如果研究對(duì)象是惰性氣體(如氬)或金屬單質(zhì)(如鐵),那么粒子(style為atomic)的屬性就是:粒子編號(hào),粒子類別,x坐標(biāo),y坐標(biāo),z坐標(biāo)。lammps中內(nèi)置了相當(dāng)多的原子類型可以查看手冊(cè)尋找自己需要的類型。這里要說(shuō)的是atom style和atom type中文翻譯都一樣,具體是啥意思。簡(jiǎn)單來(lái)說(shuō)atom style定義了大類,atom type定義了大類中的小類。比如一個(gè)水分子中,所有粒子的atom_style都是full,但是在這個(gè)style下,又區(qū)分為氧原子type和氫原子type。一個(gè)模擬中可以由多個(gè)type,但是只能有一個(gè)style。atom_style定義好后就約束了data文件中具體內(nèi)容。atom_style的語(yǔ)法是

atom_style <style> #舉例:atom_style full

粒子屬性定義完后,就可以定義粒子之間相互作用的方式了。所謂相互作用其實(shí)就是粒子之間作用力的計(jì)算公式。計(jì)算公式是通過(guò)勢(shì)能來(lái)表述的。因?yàn)閯?shì)能在空間上的梯度就是力,即f=-divu。比如重力勢(shì)能對(duì)高度求導(dǎo)就是重力。為什么要用勢(shì)能表述呢?因?yàn)閯?shì)能是標(biāo)量而力是矢量,因此使用標(biāo)量在構(gòu)建或者擬合相互作用表達(dá)式時(shí)就容易的多。而在實(shí)際計(jì)算中只需要對(duì)勢(shì)能進(jìn)行求導(dǎo)就能獲得受力。任何對(duì)象的組成粒子之間相互作用都分為非鍵結(jié)和鍵結(jié)。非鍵結(jié)就是不同通過(guò)化學(xué)鍵的相互作用,鍵結(jié)就是通過(guò)化學(xué)鍵的相互作用。有些對(duì)象只有非鍵結(jié)相互作用,如液態(tài)氬或單質(zhì)鐵。有些對(duì)象兩種都有,比如水或蛋白質(zhì)。非鍵結(jié)相互作用包括范德華和靜電。在lammps非鍵結(jié)相互作用即范德華和靜電通過(guò)pair_style和pair_coeff命令設(shè)置。pair_style和pair_coeff定義了短程相互作用,即相互作用具有截?cái)喟霃。范德華僅有短程相互作用,而靜電是既有短程又有長(zhǎng)程作用,如果體系中有靜電作用,那么定義完pair_style和pair_coeff,還要定義長(zhǎng)程靜電求解類型kspace_style。用在分子體系的鍵結(jié)相互作用分為鍵(bond),鍵角(angle),二面角(dihedral),離平面(improper)。不同的對(duì)象會(huì)包含不同的鍵結(jié)作用,比如水分子就只包括鍵和鍵角,而分子量稍大的分子會(huì)需要全部的鍵結(jié)相互作用來(lái)描述。這些相互作用構(gòu)成描述體系的力場(chǎng)。這些力場(chǎng)是從量子力學(xué)出發(fā)而得到的經(jīng)典近似。以下為定義相互作用的命令:

第七條命令設(shè)置pair_style。首先定義非鍵結(jié)相互作用。不同的力場(chǎng)有不同的pair_style,模擬中要根據(jù)對(duì)象選擇合適的pair_style。lammps集成了目前幾乎所有的pair_style。模擬時(shí)只需要根據(jù)對(duì)象進(jìn)行選擇即可。pair_style對(duì)于絕大多數(shù)模擬都需要定義。不同的pair_style具有不同的語(yǔ)法,具體可查閱手冊(cè)。pair_style中一般會(huì)定義截?cái)喟霃降拇笮。以最?jiǎn)單的包含范德華的pair_style為例,其語(yǔ)法是

pair_style lj/cut <cutoff> #舉例:pair_style lj/cut 9.0

第八條命令設(shè)置bond_style。這條命令是可選的,根據(jù)不同的對(duì)象可定。其設(shè)置簡(jiǎn)單,沒(méi)啥好說(shuō)的,就是計(jì)算兩個(gè)粒子通過(guò)化學(xué)鍵的相互作用,其語(yǔ)法為

bond_style <style> #舉例:bond_style harmonic

第九條命令設(shè)置angle_style。這條命令是可選的,根據(jù)不同的對(duì)象可定。其設(shè)置簡(jiǎn)單,沒(méi)啥好說(shuō)的,就是計(jì)算相鄰三個(gè)粒子的相互作用,其語(yǔ)法為

angle_style <style> #舉例:angle_style harmonic

第十條命令設(shè)置dihedral_style。這條命令是可選的,根據(jù)不同的對(duì)象可定。其設(shè)置簡(jiǎn)單,沒(méi)啥好說(shuō)的,就是計(jì)算四個(gè)相鄰粒子之間的相互作用,其語(yǔ)法為

dihedral_style <style> #舉例:dihedral_style opls

第十一條命令設(shè)置improper_style。這條命令是可選的,根據(jù)不同的對(duì)象可定。其設(shè)置簡(jiǎn)單,沒(méi)啥好說(shuō)的,就是計(jì)算四個(gè)相鄰粒子之間的相互作用,其語(yǔ)法為

improper_style <style> #舉例:improper_style opls

第十二條命令設(shè)置special_bonds。這條命令是可選的,當(dāng)有鍵結(jié)作用出現(xiàn)時(shí),設(shè)定有鍵結(jié)相互作用的粒子之間是否還具有pair相互作用。一般來(lái)說(shuō),是沒(méi)有的,所以這條命令一般寫(xiě)作:

special_bonds lj/coul 0.0 0.0 0.5

第十三條命令設(shè)置kspace_style。這條命令是可選的,當(dāng)有靜電作用出現(xiàn)時(shí),需要加上。不同的pair_style對(duì)應(yīng)不同的kspace_style,其語(yǔ)法:

kspace_style <style> <value> #舉例:kspace_style pppm 0.0001

以上就是第一部分的內(nèi)容。對(duì)于不同的研究對(duì)象力場(chǎng)設(shè)置部分會(huì)有所不同。如簡(jiǎn)單液體只需要pair_style就行了,分子物質(zhì)一般都會(huì)設(shè)置,晶體也一般只需要pair_style。如果你的對(duì)象有多種物質(zhì)還可以使用hybrid類的力場(chǎng)命令進(jìn)行組合設(shè)置。

第二部分 定義對(duì)象的幾何模型

第十四條命令read_data。這條命令用來(lái)讀取粒子的初始位置,鍵結(jié)信息等。這些信息寫(xiě)在一個(gè)data文件中。建立data文件稱為lammps的建模過(guò)程。lammps的建模是個(gè)非常重要的步驟。我們可以在in文件中利用creat_box和create_atoms命令建立。但是這兩條命令只能用來(lái)創(chuàng)建晶體結(jié)構(gòu)。對(duì)于分子結(jié)構(gòu)需要自己使用第三方工具或者編程建立data文件實(shí)現(xiàn)。我習(xí)慣于自己寫(xiě)data文件建立模擬對(duì)象的構(gòu)型。lammps的data文件建立比較麻煩,一般要使用多個(gè)工具聯(lián)合建模。這里不再贅述,后面會(huì)專門(mén)講lammps建立data文件的教程。data文件中的信息要根據(jù)atom_style和分子的拓?fù)浣Y(jié)構(gòu)進(jìn)行建立。以一個(gè)由分子構(gòu)成的模擬對(duì)象,也即atom_style是full為例data文件內(nèi)容為:

lammps description #第一行寫(xiě)描述信息

116803 atoms #共多少個(gè)atoms

70386 bonds #共多少個(gè)bonds

41643 angles #共多少個(gè)angles

13700 dihedrals #共多少個(gè)dihedrals

2550 impropers #共多少個(gè)impropers

191 atom types #共多少個(gè)atoms types

195 bond types #共多少個(gè)bonds types

356 angle types #共多少個(gè)angles types

548 dihedral types #共多少個(gè)dihedrals types

102 improper types #共多少個(gè)angles types

0 97.92 xlo xhi #x方向的上下邊界

0 97.92 ylo yhi #y方向的上下邊界

-15.0 160 zlo zhi #z方向的上下邊界

masses

1 15.9994 # 1號(hào)atom type的質(zhì)量

2 1.008 # 2號(hào)atom type的質(zhì)量

......

atoms

1 1 1 0.0 0.0 0.0 145.0 #粒子編號(hào),粒子所屬分子的編號(hào),粒子類別(atom type),帶電量,x坐標(biāo),y坐標(biāo),z坐標(biāo)

2 1 1 0.0 0.0 2.04 147.04

......

bonds

1 2 13825 13826 #bond的編號(hào),bond type,組成bond的第一個(gè),第二個(gè)的原子編號(hào)

2 3 13826 13827

.......

angles

1 2 13825 13826 13827 #angle的編號(hào),angle type,組成angle的第一個(gè),第二個(gè),第三個(gè)的原子編號(hào)

2 3 13825 13826 13828

......

dihedrals

1 1 13825 13826 13828 13829 #dihedral的編號(hào),dihedral type,組成dihedral的第一個(gè),第二個(gè),第三個(gè),第四個(gè)的原子編號(hào)

2 2 13826 13828 13829 13830

......

impropers

1 1 13826 13825 13827 13828 #improper的編號(hào),improper type,組成improper的第一個(gè),第二個(gè),第三個(gè),第四個(gè)的原子編號(hào)

2 2 13835 13834 13836 13837

......

第三部分 定義對(duì)象的物理模型

這部分命令用來(lái)指定勢(shì)函數(shù)中的具體參數(shù)是多少。

第十五條命令pair_coeff。根據(jù)定義的pair_style和粒子type的個(gè)數(shù)來(lái)決定pair_coeff的內(nèi)容和條數(shù)。具體書(shū)寫(xiě)時(shí)可查看手冊(cè)中對(duì)應(yīng)的pair_style的例子。如pair_style是lj/cut,粒子type個(gè)數(shù)為2,則該命令寫(xiě)為

pair_coeff 1 1 0.1553 3.166 9.0

pair_coeff 1 2 0.1423 3.125 9.0

pair_coeff 2 2 0.2753 3.196 9.0

第十六條命令bond_coeff。根據(jù)定義的bond_style和bond type的個(gè)數(shù)來(lái)決定bond_coeff的內(nèi)容和條數(shù)。具體書(shū)寫(xiě)時(shí)可查看手冊(cè)中對(duì)應(yīng)的bond_style的例子。

第十七條命令angle_coeff。根據(jù)定義的angle_style和angle type的個(gè)數(shù)來(lái)決定angle_coeff的內(nèi)容和條數(shù)。具體書(shū)寫(xiě)時(shí)可查看手冊(cè)中對(duì)應(yīng)的angle_style的例子。

第十八條命令dihedral_coeff。根據(jù)定義的dihedral_style和dihedral type的個(gè)數(shù)來(lái)決定dihedral_coeff的內(nèi)容和條數(shù)。具體書(shū)寫(xiě)時(shí)可查看手冊(cè)中對(duì)應(yīng)的dihedral_style的例子。

第十九條命令improper_coeff。根據(jù)定義的improper_style和improper type的個(gè)數(shù)來(lái)決定improper_coeff的內(nèi)容和條數(shù)。具體書(shū)寫(xiě)時(shí)可查看手冊(cè)中對(duì)應(yīng)的improper_style的例子。

第四部分 定義分組信息

很多情況下我們需要對(duì)研究對(duì)象的不同部分施加不同的操作。lammps中通過(guò)對(duì)對(duì)象中的原子進(jìn)行分組實(shí)現(xiàn)。

第二十條命令 group。group命令可以采用多種方式進(jìn)行分組定義

region fix block inf inf inf inf inf 10.0 units box

group fix region fix

group thermos type 1

group a1 id >= 76

group sub id 100:10000:10

group big molecule 123

group polya molecule <> 50 250

group boundary subtract all a2 a3

group boundary union lower upper

group boundary intersect upper flow

group boundary delete

group mine dynamic all region myregion every 100

第五部分 弛豫

在弛豫部分我們首先需要定義基本輸出日志信息,這是為了監(jiān)控我們的模擬是否正確。一般輸出的都是系統(tǒng)量,如溫度、壓強(qiáng)、總動(dòng)能、總勢(shì)能等。所以

第二十一條命令thermo。thermo命令控制每多少步輸出一次。如果不設(shè)置lammps默認(rèn)只在第一步和最后一步輸出日志信息。語(yǔ)法是

thermo <n> #舉例:thermo 1000

第二十二條命令thermo_style。thermo_style控制究竟要輸出哪些全局量。具體能夠輸出哪里量可查閱手冊(cè)thermo_style命令描述。如果不設(shè)置該條命令,那么lammps就會(huì)輸出默認(rèn)信息。其語(yǔ)法是

thermo_style <style> <args>

第二十三條命令dump。該條命令用來(lái)每隔一段時(shí)間輸出系統(tǒng)所有原子的信息,這些信息可以是原子的位置,受力,速度,質(zhì)量,編號(hào)等。當(dāng)你要輸出所有原子的信息時(shí)就是用這條命令。將這條命令的輸出文件進(jìn)行可視化可以初步判斷你的模擬是否正確。以下是最常用的設(shè)置

dump 1 all custom 10000 dump.lammpstrj id mol type mass x y z

第二十四條命令minimize。這條命令是可選的。如果你的幾何模型和物理模型都設(shè)置正確,那么可以不用這條命令。這條命令執(zhí)行的是能量最小化。為什么要執(zhí)行能量最小化呢?這是因?yàn)槲覀冊(cè)诮5臅r(shí)候,有時(shí)幾何模型和物理模型不是那么完美,經(jīng)常出現(xiàn)的就是某些原子之間的位置靠得太近,這就會(huì)導(dǎo)致系統(tǒng)的能量很大,容易出現(xiàn)運(yùn)行錯(cuò)誤。執(zhí)行能量最小化就是將系統(tǒng)的能量調(diào)節(jié)至盡可能小。lammps實(shí)現(xiàn)這個(gè)命令的方式是通過(guò)計(jì)算原子之間的受力和能量,通過(guò)一些優(yōu)化算法,移動(dòng)原子的位置從而是體系總能量達(dá)到盡可能小。所以如果你發(fā)現(xiàn)你的體系中能量超級(jí)大,壓強(qiáng)超級(jí)大而運(yùn)行出錯(cuò),那么你可以試試能量最小化命令。一般使用下面這些參數(shù)設(shè)置就能實(shí)現(xiàn)你的目的,如果能量最小化后還出錯(cuò),那么你就要排查你的幾何模型和物理模型是不是太差了,能量最小化也無(wú)能為力了。至于能量最小化采用什么算法,用默認(rèn)的就行了,也就是說(shuō)你不用有額外的設(shè)置,直接minimize就行了,花里胡哨的沒(méi)啥用。

minimize 0.0 1.0e-8 1000 100000

第二十五條命令compute。該命令是lammps中極其重要的一類命令。當(dāng)你需要計(jì)算某個(gè)量時(shí)你就可以考慮使用compute命令。lammps提供了豐富的compute命令,具體使用時(shí),根據(jù)你的需求,查閱手冊(cè)看lammps中有沒(méi)有實(shí)現(xiàn)你需要的功能的compute。compute和dump一樣可以針對(duì)某個(gè)group的粒子進(jìn)行計(jì)算。一個(gè)in文件中可以有多個(gè)compute,但每個(gè)compute的id不能重復(fù)。計(jì)算出來(lái)的值可以被很多命令使用,如dump,fix等。最常用的就是計(jì)算溫度,如

compute myt all temp

第二十六條命令velocity。velocity命令有不少功能,最常用的功能就是根據(jù)目標(biāo)溫度創(chuàng)建粒子的初始速度。最常用的形式為

velocity all create 300.0 4928459 rot yes dist gaussian

第二十七條命令timestep。該命令設(shè)定求解牛頓差分方程中時(shí)間步的大小。不同體系會(huì)要求不同大小的時(shí)間步長(zhǎng)

。一般來(lái)說(shuō)分子體系會(huì)要求時(shí)間步長(zhǎng)為1fs。簡(jiǎn)單液體體系可放寬至5fs。晶體體系可以是2fs。反應(yīng)力場(chǎng)reaxff體系會(huì)要求更小的時(shí)間步,0.1-0.5fs。時(shí)間步長(zhǎng)的選取原子是,要足夠小以準(zhǔn)確捕捉體系中最快的動(dòng)力學(xué)過(guò)程,又要足夠大以節(jié)省計(jì)算時(shí)間。通常你的時(shí)間步長(zhǎng)大小參考文獻(xiàn)設(shè)置就行了。或按照上面舉的例子進(jìn)行設(shè)置,就沒(méi)啥問(wèn)題了。

第二十八條命令fix。fix是lammps中另一類極其重要的命令。其實(shí)現(xiàn)功能就是對(duì)某個(gè)組的粒子施加某種操作。最重要的功能就是用來(lái)更新某個(gè)組的粒子的速度和位置,因此每個(gè)in文件中都必然會(huì)有fix命令。我們通常所說(shuō)的系綜就是通過(guò)fix命令實(shí)現(xiàn)。如fix nve實(shí)現(xiàn)nve系綜,fix nvt實(shí)現(xiàn)nvt系綜,fix npt實(shí)現(xiàn)npt系綜。關(guān)于系綜這里要多說(shuō)一點(diǎn)。系綜是統(tǒng)計(jì)力學(xué)的基本概念。所謂系綜,簡(jiǎn)單來(lái)說(shuō)就是系統(tǒng)無(wú)限個(gè)副本的集合。但是分子動(dòng)力學(xué)中的系綜和統(tǒng)計(jì)力學(xué)上的系綜有聯(lián)系又有區(qū)別。比如在lammps中nve系綜就是只直接求解牛頓方程而不對(duì)速度和位置做任何干涉的系綜。把fix nve與溫度控制fix聯(lián)用就實(shí)現(xiàn)nvt系統(tǒng),把fix nve與溫度控制的fix和壓強(qiáng)控制的fix同時(shí)聯(lián)用就實(shí)現(xiàn)了npt系綜。也可以只使用fix nvt命令直接實(shí)現(xiàn)nvt,只使用fix npt命令實(shí)現(xiàn)npt系綜。系綜當(dāng)你的邊界條件是p或f時(shí)系統(tǒng)體積是不變的。當(dāng)你使用s邊界時(shí)系統(tǒng)體積是可變得,此時(shí)也就沒(méi)有什么nvt或nve了。但是需要強(qiáng)調(diào)的是npt是通過(guò)改變盒子體積來(lái)實(shí)現(xiàn)控壓的所以在控壓的方向上只能使用p邊界條件。所以基本上你可以忘記系綜的概念。把fix nve或fix nvt或fix npt當(dāng)做更新原子位置和速度的命令,fix nve直接根據(jù)牛二進(jìn)行更新,fix nvt在更新的時(shí)候會(huì)對(duì)速度進(jìn)行干涉以達(dá)到設(shè)定的溫度,fix npt會(huì)對(duì)速度和位置同時(shí)進(jìn)行干涉以達(dá)到設(shè)定的溫度和壓強(qiáng)。至于nvt和nve下體積變不變?nèi)Q于你的邊界條件。同一個(gè)系統(tǒng)中可能同時(shí)存在fix nve和fix nvt,這是因?yàn)橛行┠M要求一部分原子的位置和速度直接根據(jù)牛二進(jìn)行更新,而另一部分卻要求在更新原子速度時(shí)保持溫度不變。但是不會(huì)同時(shí)出現(xiàn)fix nve和fix npt,因?yàn)閒ix npt在根據(jù)設(shè)定的壓強(qiáng)在更新原子位置和速度時(shí),計(jì)算的是所有原子對(duì)壓強(qiáng)的貢獻(xiàn)。所以fix npt不適合局部控壓。

某個(gè)組的粒子只有使用了fix+系綜的命令才會(huì)更新位置和速度。因此如果你想保持某組的粒子在模擬過(guò)程中固定不動(dòng),最簡(jiǎn)單的辦法就是不對(duì)它施加任何fix命令。通常我們?cè)谑褂胒ix時(shí)都會(huì)使用先來(lái)一條fix+系綜,然后調(diào)用其他fix實(shí)現(xiàn)某些目的,如要讓模擬對(duì)象發(fā)生一定變形,就用fix deform,要讓對(duì)對(duì)象施加一個(gè)電場(chǎng)就用fix efield。當(dāng)你想用實(shí)現(xiàn)某種操作時(shí),打開(kāi)手冊(cè)查看fix的一些命令看有沒(méi)有實(shí)現(xiàn)你功能的fix命令。fix的某些命令也可以用于輸出信息,如fix ave/chunk輸出模擬盒子不同位置處的原子屬性的時(shí)間平均值;fix ave/time輸出某些全局量的平均值。

第二十九條命令 run。終于要開(kāi)始跑了。這個(gè)命令最常用的就是你要run多少步。在弛豫階段一般會(huì)run 1ns的長(zhǎng)度,通常是1000000步。這里可以控制采用什么算法對(duì)牛二進(jìn)行求解。采用默認(rèn)的就行,也就是說(shuō)你不用額外設(shè)置,直接run就行了;ɡ锖诘臎](méi)啥用,默認(rèn)的往往是最好的。其設(shè)置是

run 1000000

以上就是弛豫部分。為什么要進(jìn)行弛豫呢?這是因?yàn)橄到y(tǒng)中原子的初始位置和初始速度都是認(rèn)為設(shè)定的,這是系統(tǒng)的狀態(tài)并不在熱力學(xué)平衡狀態(tài)。通常我們希望在正式模擬時(shí)系統(tǒng)的起點(diǎn)在熱力學(xué)平衡狀態(tài)。只有經(jīng)過(guò)弛豫后系統(tǒng)才會(huì)達(dá)到熱力學(xué)平衡狀態(tài)。在弛豫開(kāi)始的時(shí)候你會(huì)看到系統(tǒng)的溫度和壓強(qiáng)會(huì)有非常劇烈的變化,這是正常的。經(jīng)過(guò)一段時(shí)間后溫度和壓強(qiáng)這些量就會(huì)在一個(gè)值附近小幅振蕩。那么怎么判斷系統(tǒng)是否達(dá)到平衡了呢?這是個(gè)挺麻煩的事,不同的體系有不同的標(biāo)準(zhǔn)。通常我們會(huì)看系統(tǒng)的溫度,壓強(qiáng),能量這些值不會(huì)在持續(xù)上升或下降就差不多平衡了。但是這個(gè)事很難說(shuō)。所以為了保險(xiǎn)起見(jiàn),通常模擬中會(huì)設(shè)置冗余的弛豫步數(shù)來(lái)確保系統(tǒng)達(dá)到平衡。這個(gè)值一般是1000000步。需要說(shuō)明的是,系統(tǒng)溫度是根據(jù)系統(tǒng)動(dòng)能算出來(lái)的,一般很快就不咋變了。但是壓強(qiáng)是根據(jù)受力算出來(lái)的振蕩會(huì)很大。一般要看溫度或壓強(qiáng)的時(shí)間平均值,也就是不同時(shí)間步的平均值來(lái)判斷系統(tǒng)的溫度或壓強(qiáng)是否達(dá)到了設(shè)定值。

第三十條命令write_restart。這個(gè)命令可以寫(xiě)出restart文件,以便后續(xù)從這個(gè)點(diǎn)接著算。restart文件是一個(gè)二進(jìn)制文件。使用這條命令是個(gè)好的習(xí)慣,可以讓你省很多事,又不會(huì)額外增加計(jì)算量。其語(yǔ)法是

write_restart <filename> #舉例write_restart restart.equ

最后一部分 正式模擬

這部分其實(shí)才是模擬的核心。不同的研究對(duì)象和目的設(shè)置和使用的命令都會(huì)不同。無(wú)法在基本流程中寫(xiě),將在后面的專題中介紹。

5. 建模和后處理軟件概述

本節(jié)介紹了包裝 lammps 的商業(yè)和免費(fèi)軟件,為開(kāi)發(fā)模型、運(yùn)行仿真和分析結(jié)果提供了一個(gè)用戶友好的環(huán)境。

materials design, inc.

materials design, inc. 開(kāi)發(fā)了medea®,這是一種原子模擬和建模環(huán)境,可提供與lammps一起使用的生產(chǎn)力、模型構(gòu)建和分析工具。medea®使用流程圖簡(jiǎn)化了lammps模擬。生成的lammps模擬很容易在同事之間共享、編輯以供將來(lái)重復(fù)使用,并且可以由lammps專家進(jìn)行定制。

使用medea®,可以使用各種晶體系統(tǒng)方法構(gòu)建原子模型,并使用medea®非晶和熱固性構(gòu)建器,這是為lammps創(chuàng)建輸入的工具。它們還為預(yù)測(cè)機(jī)械性能提供了經(jīng)過(guò)驗(yàn)證的方法,包括基于lammps模擬的彈性常數(shù)、擴(kuò)散率、傳輸性能和內(nèi)聚能。medea®還有助于管理用于模擬有機(jī)、無(wú)機(jī)和金屬系統(tǒng)的力場(chǎng);材料設(shè)計(jì)團(tuán)隊(duì)在為有機(jī)和金屬系統(tǒng)開(kāi)發(fā)精確力場(chǎng)方面擁有良好的記錄。

scienomics

scienomics公司開(kāi)發(fā)了一個(gè)lammps接口,作為其材料和工藝模擬 (maps) 平臺(tái)的一部分,該平臺(tái)允許lammps新手和專家快速有效地創(chuàng)建lammps輸入文件-用于原子模擬和dpd模擬。scienomics為 lammps用戶提供完整的電話和電子郵件支持。maps中的lammps插件還允許用戶為lammps創(chuàng)建輸入文件,對(duì)lammps模擬的輸出進(jìn)行可視化和分析。maps也有許多力場(chǎng)可供選擇(用于原子建模和粗粒度建模)與 lammps。他們還有一個(gè)名為amorphous builder和crosslink builder的工具,可用于為lammps創(chuàng)建輸入幾何圖形,以及用于分析reaxff模擬結(jié)果的工具和其他強(qiáng)大的功能,用于計(jì)算彈性屬性、結(jié)構(gòu)屬性、傳輸屬性等關(guān)鍵屬性.

maps是一個(gè)綜合性的軟件平臺(tái),結(jié)合了構(gòu)建-模擬-分析的工作流程。這使用戶能夠以高效和無(wú)縫的方式對(duì)復(fù)雜材料系統(tǒng)進(jìn)行建模,從包括lammps在內(nèi)的一系列仿真引擎中進(jìn)行選擇以對(duì)過(guò)程進(jìn)行建模,并分析輸出以提取許多屬性,例如機(jī)械、熱、傳輸、流變等. 有關(guān)所有這些工具和服務(wù)的更多詳細(xì)信息,請(qǐng)參見(jiàn)scienomics網(wǎng)頁(yè)。

xenoview

由sergei shenogin倫斯勒納米技術(shù)中心)開(kāi)發(fā)。該軟件對(duì)非商業(yè)用戶免費(fèi)。 xenoview是基于windows的分子動(dòng)力學(xué)模擬軟件。它的界面為從頭開(kāi)始構(gòu)建結(jié)構(gòu)提供了廣泛的工具。此外,您可以從多種格式(例如pdb格式)導(dǎo)入結(jié)構(gòu)。 鍵可以自動(dòng)定義,力場(chǎng)可以自動(dòng)分配。xenoview可以導(dǎo)出可在lammps中用于大規(guī)模模擬的數(shù)據(jù)文件和輸入腳本。有關(guān)詳細(xì)信息,請(qǐng)參閱xenoview主頁(yè)。

atomman

atomman 工具是由nist的lucas hale開(kāi)發(fā)的開(kāi)源(免費(fèi))軟件。它是一個(gè)用于與大規(guī)模原子系統(tǒng)交互的python包。它允許用戶完全從python中準(zhǔn)備、運(yùn)行和分析md模擬。

l 允許對(duì)數(shù)百萬(wàn)個(gè)原子進(jìn)行高效和快速的計(jì)算,每個(gè)原子都有許多自由定義的每個(gè)原子屬性。

l 生成包含缺陷的原子系統(tǒng),例如點(diǎn)缺陷和位錯(cuò)單極子。

l 用于缺陷分析的內(nèi)置函數(shù),例如滑移向量和 nye 張量。

l 直接從 python 調(diào)用 lammps 并立即檢索結(jié)果數(shù)據(jù)或 lammps 錯(cuò)誤語(yǔ)句。

l 輕松將系統(tǒng)與 ase 和 pymatgen 的其他 python 原子環(huán)境相互轉(zhuǎn)換。

l 可以加載基于 cif 晶體結(jié)構(gòu)文件、poscar 和 lammps 原子和轉(zhuǎn)儲(chǔ)文件的系統(tǒng)。

l 內(nèi)置單位轉(zhuǎn)換工具。

l 獨(dú)立于平臺(tái):適用于 linux、windows 和 mac。

loos = lightweight object-oriented structure analysis library

loos是一個(gè)用于分析分子動(dòng)力學(xué)模擬的軟件包。 它與包無(wú)關(guān)、開(kāi)源,可在所有主要的linux發(fā)行版以及osx上運(yùn)行。它分發(fā)了大約150個(gè)預(yù)打包的分析工具,范圍從標(biāo)準(zhǔn)任務(wù)(軌跡操作、主成分分析等)到新工具(評(píng)估模擬收斂、測(cè)量膜和膜蛋白的特性)。 此外,它專為快速開(kāi)發(fā)新的分析工具而設(shè)計(jì),尤其是在使用python包裝器時(shí)。loos可從https://github.com/grossfieldlab/loos下載

freud = analysis library for post-processing md and mc data

freud是一個(gè)python庫(kù),由sharon glotzer的小組(密歇根大學(xué))開(kāi)發(fā)和支持,它提供了一套簡(jiǎn)單、靈活、強(qiáng)大的工具,用于分析從分子動(dòng)力學(xué)或蒙特卡羅模擬中獲得的軌跡。 高性能、并行化 c++ 用于計(jì)算標(biāo)準(zhǔn)工具,例如徑向分布函數(shù)、相關(guān)函數(shù)、階參數(shù)和聚類,以及包括平均力和扭矩 (pmft)勢(shì)和局部環(huán)境匹配在內(nèi)的原始分析方法。freud庫(kù)支持多種輸入格式和輸出numpy數(shù)組,支持與許多典型材料科學(xué)工作流程的科學(xué)python生態(tài)系統(tǒng)集成。

教程:g(r) 分析

文檔:https://freud.readthedocs.io

源代碼:https://github.com/glotzerlab/freud

molecular builders

要模擬分子系統(tǒng),lammps要求您輸入分子拓?fù)浣Y(jié)構(gòu)(鍵、角、二面體等列表)以及適合您的模型的力場(chǎng)系數(shù)。因此,構(gòu)建分子系統(tǒng)的任務(wù)是一個(gè)預(yù)處理步驟,并且本身可能是一項(xiàng)復(fù)雜的任務(wù)。

本節(jié)介紹有助于自動(dòng)化此過(guò)程的工具。

一般來(lái)說(shuō),封裝可以從鍵拓?fù)渫茢嘟嵌、二面角和不正確的相互作用。 他們有生成分子幾何的命令。他們可以從packmol和其他pdb文件生成器生成的文件中讀取坐標(biāo)。

enhanced monte carlo (emc)

由位于https://verizon.net的't veld (basf)的pieter j.開(kāi)發(fā)和維護(hù)。

增強(qiáng)型monte carlo或簡(jiǎn)稱emc為使用compass、charmm、opls、martini、dpd或膠體力場(chǎng)創(chuàng)建和操作粒子模擬的輸入結(jié)構(gòu)提供了環(huán)境。為此,腳本語(yǔ)言管理對(duì)其功能的訪問(wèn)。當(dāng)前版本通過(guò)smiles字符串提供對(duì)分子或粗粒度結(jié)構(gòu)的操作,在需要時(shí)為選定的力場(chǎng)鍵入這些結(jié)構(gòu),并應(yīng)用蒙特卡羅原理構(gòu)建構(gòu)象以不重疊原子。emc提供lammps、pdb和xyz格式的輸出端口?梢栽https://montecarlo.sourceforge.net上找到適用于linux、macos或windows的編譯版本。

moltemplate

由andrew jewett (ucsb) 開(kāi)發(fā)和維護(hù)。

該工具與 lammps 一起在tools/moltemplate目錄中分發(fā)。有關(guān)更多詳細(xì)信息,請(qǐng)參見(jiàn)該目錄和moltemplate主頁(yè)。

moltemplate專為構(gòu)建粗粒度生物分子模型而設(shè)計(jì)。moltemplate可以同時(shí)創(chuàng)建:lammps data文件(包含幾何和拓?fù)洌┖蚻ammps input腳本(包含力場(chǎng)、修復(fù)和組)。與其他轉(zhuǎn)換工具生成的文件不同,moltemplate允許用戶訪問(wèn) lammps中可用的所有力場(chǎng)。用戶可以將分子保存為moltemplate緊湊、可讀的模板文件格式(“.lt”),并與他人共享。分子可以用作更大分子的構(gòu)建塊。 “罐裝”力場(chǎng)(例如 dreiding、gaff、trappe 和用戶創(chuàng)建)也可以(原則上)以這種格式保存并稍后應(yīng)用于分子。

分子可以被復(fù)制、組合和鏈接在一起以定義新分子。(這些可用于定義更大的分子。支持無(wú)限級(jí)別的對(duì)象組合、嵌套和繼承。)一旦構(gòu)建,可以自定義單個(gè)分子和子單元(原子和鍵,以及子單元可以移動(dòng)、刪除和替換)。

vmd topotools

由axel kohlmeyer(temple u開(kāi)發(fā)和維護(hù)。

有關(guān)詳細(xì)信息,請(qǐng)參閱 topotools 主頁(yè)。

topotools是一個(gè)分子生成器,它利用vmd和tcl的強(qiáng)大功能創(chuàng)建lammps data文件并將它們轉(zhuǎn)換為其他格式或從其他格式轉(zhuǎn)換。topotools有兩個(gè)組件:一個(gè)可以提取和操作拓?fù)湫畔⒌闹虚g件腳本,以及在它之上構(gòu)建的幾個(gè)高級(jí)應(yīng)用程序,例如可以使其能夠讀取/寫(xiě)入數(shù)據(jù)文件、復(fù)制和合并系統(tǒng)。 與vmd一起,topotools可以從pdb文件、psf文件和原子對(duì)距離推斷拓?fù)浣Y(jié)構(gòu),使蛋白質(zhì)溶劑化。

intermol a conversion tool for molecular dynamics simulations

intelmol是用python編寫(xiě)的,可以在intermol中本地進(jìn)行 desmond<=>gromacs<=>lammps轉(zhuǎn)換;amber->x 是通過(guò)將amber 轉(zhuǎn)換為 gromacs,然后使用parmed 轉(zhuǎn)換為其他程序來(lái)執(zhí)行的。amber->charmm 由parmed直接執(zhí)行。

avogadro

作為開(kāi)源項(xiàng)目開(kāi)發(fā)和維護(hù)。有關(guān)詳細(xì)信息,請(qǐng)參閱 avogadro 主頁(yè)。

avogadro是一款先進(jìn)的分子編輯器和可視化工具,專為計(jì)算化學(xué)、分子建模、生物信息學(xué)、材料科學(xué)和相關(guān)領(lǐng)域的跨平臺(tái)使用而設(shè)計(jì)。 它提供靈活的高質(zhì)量渲染和強(qiáng)大的插件架構(gòu)。

packmol

作為開(kāi)源項(xiàng)目開(kāi)發(fā)和維護(hù)。有關(guān)詳細(xì)信息,請(qǐng)參閱 packmol 主頁(yè)。

packmol通過(guò)在定義的空間區(qū)域中包裝分子來(lái)創(chuàng)建分子動(dòng)力學(xué)模擬的初始點(diǎn)。包裝保證短程排斥相互作用不會(huì)破壞模擬。

atomsk

作為開(kāi)源項(xiàng)目開(kāi)發(fā)和維護(hù)。有關(guān)詳細(xì)信息,請(qǐng)參閱 atomsk 主頁(yè)。

atomsk旨在創(chuàng)建、操作和轉(zhuǎn)換原子系統(tǒng)。它支持多種文件格式,其中包括 lammps文件格式,還支持vasp、quantum espresso、imd、dl_poly、atomeye cfg格式或xcrysden xsf格式,可以輕松轉(zhuǎn)換文件以進(jìn)行從頭計(jì)算、經(jīng)典勢(shì)模擬或可視化。 此外,atomsk還可以對(duì)原子位置進(jìn)行一些簡(jiǎn)單的變換,如旋轉(zhuǎn)、變形、插入位錯(cuò)。

octa 和j-octa

octa是一個(gè)開(kāi)源軟件包,由模擬引擎(分子動(dòng)力學(xué)、流變模擬、自洽場(chǎng)理論、有限元方法等)和用于建模軟物質(zhì)系統(tǒng)的gui(可視化、簡(jiǎn)單分子構(gòu)建器和分析工具)組成。

octa還為多種模擬器的協(xié)作使用提供了環(huán)境,即多物理場(chǎng)和多尺度模擬。

lammps和cognac(md 引擎)文件之間的轉(zhuǎn)換程序也包括在內(nèi)。此功能使lammps用戶能夠使用octa gui并結(jié)合其他理論(如 scft)運(yùn)行經(jīng)典的md模擬。

通過(guò)使用商業(yè)版本的j-octa,還可以進(jìn)行全原子和粗粒度md的復(fù)雜分子構(gòu)建。

molecular simulation design framework (mosdef)

mosdef提供了一組可擴(kuò)展的python工具,旨在促進(jìn)使用分子模擬對(duì)軟物質(zhì)系統(tǒng)進(jìn)行初始化、原子分型和篩選。

automated topology builder (atb) and repository - molecular topology building blocks for organic molecules

atb項(xiàng)目由alan e. mark教授(澳大利亞昆士蘭大學(xué))領(lǐng)導(dǎo);atb網(wǎng)站上有一份其他貢獻(xiàn)者的列表。

atb以與lammps和其他分子動(dòng)力學(xué)軟件包兼容的格式提供有機(jī)分子的拓?fù)湮募ammps拓?fù)湮募С质褂胢oltemplate工具(隨lammps分發(fā))構(gòu)建復(fù)雜系統(tǒng)。moltemplate允許結(jié)合力場(chǎng)參數(shù)和分子模板文件來(lái)構(gòu)建復(fù)雜的系統(tǒng),從而實(shí)現(xiàn)類似于在生物分子md模擬包(例如gromos、amber、charmm等)中構(gòu)建系統(tǒng)拓?fù)涞墓ぷ髁鞒獭?br />
atb 站點(diǎn)為各種分子(> 20,000并且還在增長(zhǎng))提供拓?fù)浣Y(jié)構(gòu)。可以使用基于名稱或化學(xué)式的內(nèi)部搜索工具并從結(jié)果列表中選擇分子來(lái)找到分子。一旦進(jìn)入感興趣分子的頁(yè)面,您可以選擇“molecular dynamics (md) files”并選擇“l(fā)ammps”作為輸出“格式”如果數(shù)據(jù)庫(kù)中不存在分子,可以通過(guò)提交新的要處理的分子。

data sites

本節(jié)列出了提供可與lammps一起使用的數(shù)據(jù)的www站點(diǎn)。例如,原子配置,可用于初始化lammps模擬的模型。或者使用lammps完成的計(jì)算存檔,可以瀏覽以獲取物理或模型。

jarvis database for md potential calculations on dft geometries

jarvis for force-fields (jarvis-ff)是一個(gè)高通量計(jì)算數(shù)據(jù)庫(kù),用于對(duì)具有各種力場(chǎng)/原子間勢(shì)的密度泛函理論 (dft)優(yōu)化的幾何結(jié)構(gòu)進(jìn)行l(wèi)ammps計(jì)算。 該項(xiàng)目的目標(biāo)是通過(guò)網(wǎng)絡(luò)界面為評(píng)估力場(chǎng)提供一個(gè)簡(jiǎn)單的查找表,并提高數(shù)據(jù)的可重復(fù)性。jarvis-ff是美國(guó)國(guó)家標(biāo)準(zhǔn)與技術(shù)研究院 (nist)材料基因組計(jì)劃(mgi)的一部分。

orsi group at queen mary university of london

mario orsi的小組為各種水系統(tǒng)(不同的力場(chǎng))收集了很好的lammps輸入腳本和數(shù)據(jù)文件。他們也在為蛋白質(zhì)和脂質(zhì)膜系統(tǒng)添加文件。

[ Last edited by 月只藍(lán) on 2022-5-22 at 17:11 ]
回復(fù)此樓

» 收錄本帖的淘帖專輯推薦

缺陷計(jì)算 科研軟件 分子動(dòng)力學(xué)模擬相關(guān) lammps

» 本帖已獲得的紅花(最新10朵)

» 猜你喜歡

» 本主題相關(guān)商家推薦: (我也要在這里推廣)

» 搶金幣啦!回帖就可以得到:

查看全部散金貼

已閱   回復(fù)此樓   關(guān)注TA 給TA發(fā)消息 送TA紅花 TA的回帖

學(xué)員0C9LnR

鐵蟲(chóng) (小有名氣)



luyaobao(金幣+5): 謝謝參與
點(diǎn)贊!期待專題六!
41樓2022-11-21 15:42:16
已閱   回復(fù)此樓   關(guān)注TA 給TA發(fā)消息 送TA紅花 TA的回帖
查看全部 64 個(gè)回答

luwis

至尊木蟲(chóng) (職業(yè)作家)



luyaobao(金幣+5): 謝謝參與
2樓2022-05-22 15:27:52
已閱   回復(fù)此樓   關(guān)注TA 給TA發(fā)消息 送TA紅花 TA的回帖

whgregory

木蟲(chóng) (正式寫(xiě)手)



luyaobao(金幣+5): 謝謝參與
5樓2022-05-22 17:16:48
已閱   回復(fù)此樓   關(guān)注TA 給TA發(fā)消息 送TA紅花 TA的回帖
簡(jiǎn)單回復(fù)
tzynew6樓
2022-05-22 18:37   回復(fù)  
luyaobao(金幣+5): 謝謝參與
i 發(fā)自小木蟲(chóng)Android客戶端
2022-05-22 19:31   回復(fù)  
luyaobao(金幣+5): 謝謝參與
yes 發(fā)自小木蟲(chóng)Android客戶端
jiaoxg12樓
2022-05-23 09:14   回復(fù)  
luyaobao(金幣+5): 謝謝參與
2022-05-23 17:13   回復(fù)  
luyaobao(金幣+5): 謝謝參與
5 發(fā)自小木蟲(chóng)IOS客戶端
提示: 如果您在30分鐘內(nèi)回復(fù)過(guò)其他散金貼,則可能無(wú)法領(lǐng)取此貼金幣
最具人氣熱帖推薦 [查看全部] 作者 回/看 最后發(fā)表
[考研] 085700資環(huán)求調(diào)劑,初始279,六級(jí)已過(guò),英語(yǔ)能力強(qiáng) +4 085700資環(huán)調(diào)劑 2026-03-03 5/250 2026-03-08 13:35 by 30660438
[考研] 材料調(diào)劑 +7 ounce. 2026-03-05 15/750 2026-03-08 09:39 by wangjihu
[考研] 調(diào)劑的同學(xué),走過(guò)路過(guò),不要錯(cuò)過(guò) +8 likeihood 2026-03-06 13/650 2026-03-08 09:35 by likeihood
[考研] 【求調(diào)劑】293分環(huán)境工程求調(diào)劑材料/化工,服從調(diào)劑,抗壓能力強(qiáng)! +10 xiiiia 2026-03-04 10/500 2026-03-08 06:38 by 劉兵
[考研] 材料科學(xué)(0805)338 求調(diào)劑 +7 xiaokang3286 2026-03-07 7/350 2026-03-08 03:42 by lfhuang
[考研] 一志愿鄭大071000分?jǐn)?shù)282求調(diào)劑 +3 研研顏 2026-03-05 7/350 2026-03-07 22:33 by 帆船哥
[考研] 求調(diào)劑,不管什么專業(yè),我是可塑造的人才一枚,希望遇到知己老師撈撈我 +4 13102137290 2026-03-06 5/250 2026-03-07 21:21 by lissomchan
[考研] 085600材料與化工 298 調(diào)劑 +11 小西笑嘻嘻 2026-03-03 11/550 2026-03-07 10:30 by 一直走不要停
[考研] 求調(diào)劑 一志愿蘇州大學(xué),0856化工323分 | 本科應(yīng)化 | 有專利/競(jìng)賽/科研助手經(jīng)歷 | +5 橙子cyx 2026-03-06 6/300 2026-03-06 23:13 by L135790
[考研] 306求調(diào)劑 +7 Bahati 2026-03-05 7/350 2026-03-06 22:11 by 星空星月
[考研] 材料277分求調(diào)劑 +13 飯飯星球 2026-03-04 14/700 2026-03-06 16:10 by @颯颯颯颯
[考研] 沒(méi)上岸的看過(guò)來(lái) +3 tangxiaotian 2026-03-01 5/250 2026-03-05 13:56 by tangxiaotian
[考研] 0856材料專碩274能調(diào)劑去哪里? +3 22735 2026-03-04 4/200 2026-03-05 09:06 by 斬魂滴兔子!
[考研] 085601 材料工程 320 +6 和樂(lè)瑤 2026-03-03 6/300 2026-03-04 16:01 by chixmc
[考研] 本科太原理工采礦工程,求調(diào)劑 +3 onlx 2026-03-01 3/150 2026-03-04 15:57 by Stephen_ym
[考研] 0854總分272 +5 打小就是老實(shí)人 2026-03-02 6/300 2026-03-04 01:41 by ouhaiyu
[考研] 26考研報(bào)考西工大材料308分求調(diào)劑 +4 weizhong123 2026-03-01 5/250 2026-03-03 12:22 by weizhong123
[論文投稿] 通訊作者寫(xiě)誰(shuí),問(wèn)題是你意想不到的問(wèn)題 15+3 阿爾法啊 2026-03-01 3/150 2026-03-03 09:13 by 北京萊茵潤(rùn)色
[考研] 261求調(diào)劑 +3 陸lh 2026-03-01 3/150 2026-03-02 19:32 by zhukairuo
[考研] 275求調(diào)劑 +3 L-xin? 2026-03-01 6/300 2026-03-02 10:22 by 熱情沙漠
信息提示
請(qǐng)?zhí)钐幚硪庖?jiàn)
天天看片天天摸天天操| 久久av色噜噜ai换脸| 情趣视频在线观看91| 亚洲第一中文字幕成人| 天天操天天干加勒比久久| 妈妈的朋友中字在线免费观看| 天天透天天舔天天操| 精久久久久久久久久久久| 夫亡人妻被强干中文字幕| 亚洲午夜高清在线观看| 丰满少妇_区二区三区| 欧美大胆a级视频秒播| 亚洲av综合av一去二区三区| 97精品国产91久久久| 青青草一个释放的网站| 亚洲精品一区二区gif| 免费24小时人妻视频| 99 re国产精品| av 一区二区三区 熟女| 中文在线字幕免费观看日韩视频| 久久热在线免费观看| 在线观看免费啪啪啪| 日韩欧美一区二区三区免费看| 日本五六十路熟女视频| 免费在线观看视频啪啪| 极品少妇高潮喷水日出白浆| 18岁禁一二三区免费体验| 黑人大吊大战亚洲女人。| 蜜臀一区二区日韩美女少妇视频| 九九热视频1这里只有精品| 久久热在线免费观看| 精品精品精品精品精品污污污污| 天天弄天天草天天日天天| 日日夜夜免费视频精品| 欧美日本亚欧在线观看| 日韩少妇免费在线播放| 亚洲综合天堂av网站在线观看| 乌克兰美女操逼高清内射视频| av男人站在线观看| 荣立三等功退休有什么待遇| 男女69视频在线观看免费| 亚洲高清一区二区三区久久| 三级欧美日韩一区二区三区| 68福利精品在线视频| 欧美aaaa性bbbbaaaa| 亚洲制服丝袜网站中文字幕| 正在播放麻豆精品一区二区| 在线 激情 亚洲 视频| 午夜五十路久久福利| 亚洲美女a级黄色在线播放| 成人资源中文在线观看| 鸡巴在里面福利视频在线观看| 人妻色综合aaaaaa网| 4438x亚洲最大的成人| 亚洲一区二区精品在线播放| 亚洲av三级电影在线观看| 久久99嫩草99久久精品| 一区二区欧美 国产日韩| 精产国品一二三产品区别97| av网页免费在线观看| 午夜福利国产精品久久久久| 中文字幕 首页 人妻| 亚洲成年人精品国产| 精品国产人伦一区二区三区| 麻豆国产精品777777在| 人妻激情综合久久久久蜜桃 | 成人免费电影二区三区 | 91精品一区一区三区| 中文字幕人妻一区二区视频系列 | 久久久亚洲熟女一区二区| 亚洲人人爽人人澡起碰av| av中文字幕国产精品| 精品人妻在线激情视频| 北野中文字幕一区二区| 69av精品国产探花| 91精品国产成人久久久久久| 日韩一区二区在线播放观看| 天天操,天天射,天天爽| 日本少妇人妻中文在线| 欧美精品一区二区三区观看| 国产熟妇色xxⅹ交白浆视频| 亚洲欧美综合另类最新| 成人人妻h在线观看| 91亚洲最新蜜桃在线| 中文字幕观看中文字幕免费| 亚洲欧洲一区二区三区在线| 网友自拍第一页99热| 九九九九九久久久国产| 9999久久久久老熟妇二区| 日本韩国欧美在线视频| 精品精品精品精品精品污污污污| 久久久久久久精品乱码| 亚洲精品9999蜜桃| 人妻色综合aaaaaa网| 69av精品国产探花| 91porny九色视频偷拍| 高潮喷水在线视频观看| 美国十次了亚洲天堂网国产| 天天日天天干天天日天天干天天 | 1级黄色片在线观看| 五月天色婷婷狠狠爱| 四虎精品久久免费最新| 天天操天天射天天操天天日| 人妻人妻在线视频网站| 亚洲春色av中文字幕| 亚洲欧美日韩中文在线观看| 不卡高清一区二区三区| 丰满人妻熟女aⅴ一区| 漂亮人妻口爆久久精品| 亚洲黄色成人一级片| 国产天堂av不卡网| 日韩人妻中文字幕二区| 亭亭五月天在线观看| 夜夜人人干人人爱人人操| 欧美日本亚欧在线观看| 北野中文字幕一区二区| 最近最新欧美日韩精品| 亚洲18片综合国产av| 亚洲另类欧美综合久久| 啪啪啪网站免费在线看| 青青青免费手机视频在线观看| 日韩成人精品久久久免费看| 琪琪日本福利伦理视频| 亚洲中文字幕最新地址| 免费24小时人妻视频| 天天想要天天操天天干| 国产资源在线观看二区| 日本一道中文字幕99| 网站在线观看蜜臀91| 视频免费在线观看网站| 猫咪亚洲中文在线中文字幕| 女人的天堂 av在线| 亚洲第一中文字幕成人| 欧美在线视频不卡一区| 熟妇人妻av无码中文字幕| 欧美肥妇久久久久久| 亚洲一区亚洲二区成人福利| 视频自拍偷拍视频自拍| 美女福利网站在线播放| 在线观看视频免费一区二区三区| 日韩激情文学在线视频| 色就色综合偷拍区欧美在线| 亚洲国产精品青青草| 亚洲欧美激情国产综合久久久| 大尺度久久久久久久| 高清av在线婷一区二区色日韩| 日本高清有码在线视频| 中国特黄色性生活片| 亚洲少妇视频在线观看| 五十岁熟女高潮喷水| 全彩漫画口工18禁| 亚洲成人五月婷婷久久综合| 午夜宅男电影av网站| 亚洲av毛片在在线播放| 亚洲中文字幕在线视频观看二区| 日韩av熟妇在线观看| 国产精品黄色片大全| 国产福利小视频在线观看网站| 黑鸡巴肏少妇逼视频| 91精品资源在线观看| 99精品视频在线在线观看| 亚洲欧美激情久久久| 亚洲精品综合欧美精品综合| 天天干夜夜撸天天操| 伊人综合在线视频免费观看| 亚洲唯美激情综合四射| 特级aaaaa黄色片| 在线 激情 亚洲 视频| 91精品国产成人久久久久久| 黄色av网址在线播放| 99999久久久精品| 日本黄页在线观看视频| 黄在线看片免费人成视频| 亚洲va999天堂va| 大乳人妻一区二区三区| 国产精品igao为爱寻找激情| 青青青免费手机视频在线观看| 天天操天天干加勒比久久| 两个人在一起靠逼啊啊啊| 国产成人情侣av在线| 92麻豆一区二区三区| 在线有码人妻自拍视频| 91日本精产品一区二区三区 | 亚洲精品9999蜜桃| 在线播放 日韩 av| 99色在线观看免费观看| 欧美猛少妇色ⅹⅹⅹⅹⅹ猛叫| 欧美成人性生活视频播放| 伊人精品久久一区二区| 老熟女xxxⅹhd老熟女性| 欧美在线观看一区二区不卡| 欧美日韩一区二区三区成人影院| 一区二区三区 国产日韩欧美| 最新免费在线观看污视频| 玖玖资源站在线观看亚洲| 日本欧美亚洲国产啊啊啊| 国产极品气质外围av| 国产又粗又长又大视频| 大尺度av毛片在线网址| 乌克兰美女操逼高清内射视频| 91香蕉国产亚洲一二三区| 亚洲综合第一区二区| 七色福利视频在线观看| 一区二区欧美 国产日韩| 顶级欧美色妇4khd| 午夜五十路久久福利| 中文字幕在线观看av观看| 亚洲第一成年偷拍视频| av资源中文字幕在线观看 | av在线观看视频免费| av无限看熟女人妻另类av| avjpm亚洲伊人久久| 国产igao激情在线视频入口| 蜜乳av一区二区三区免费观看| 亚洲精品综合欧美精品综合| 91偷拍被偷拍在线播放| 天天弄天天草天天日天天| 国产高清视频www夜色资源| 女人的天堂av在线网| 久久人妻诱惑我视频| 一区二区三区资源视频| 美女一区二区四区六区八区| 日韩久久九九精品视频| jizzjizz国产精品传媒| 国产一级一国产一级毛片 | 亚洲午夜熟女在线观看| 68福利精品在线视频| 蜜乳视频一区二区三区| 日韩精品欧美一区二区| 九九热视频1这里只有精品| 日本国产亚洲欧美色综合| av一区二区三区蜜桃| 91九色人妻在线播放| 久久无码高清免费视频| 日日躁夜夜躁狠狠操| iga肾三级算严重吗| 人妻免费视频黄片在线视频| 欧美亚洲愉拍一区二区三区| 91性高湖久久久久久久久久| 国产精品久久久99| 午夜精品久久久久久久精品乱码 | 川上优所有中文字幕在线| 老司机伊人99久久精品| 人妻激情综合久久久久蜜桃 | 国产精品久久人人添| 一区二区三区资源视频| 快色视频在线观看免费| 91九色pony蝌蚪| 欧美日本在线免费视频| 蜜桃臀av在线一区二区| 亚洲乱码av一区二区蜜桃av| 中文字幕熟女人妻一区| av在线观看视频免费| 中国精品人妻一区二区| 亚洲国产综合久久精品| 欧美大鸡吧男操女啊啊啊视频| 熟妇人妻丰满久久久久久久| 超碰在线免费观看视频97| 中文字幕一区二区三区久久久| 网友自拍第一页99热| 国产漂亮白嫩美女在线图片 | 人妻少妇视频系列视频在线| 天天做天天日天天搞| 加勒比不卡在线视频| 果冻麻豆一区二区三区| 人人妻人人爽人人爽欧美一区| 国内自拍第一区二区三区| 国产一级一国产一级毛片| 成年男女免费视频网站无毒| 男人和女人的逼视频| 国产夫妻视频在线观看免费| 91青青青国产免费高清| 欧美国产精品久久久免费| 自拍偷自拍亚洲精品10p| 18禁网站在线点击观看| 久久久久久a女人处女| 制服丝袜 中文字幕 日韩| 一区二区欧美 国产日韩| 全球高清中文字幕av| 91九色pony蝌蚪| 操死你美女在线视频| 亚洲最大先锋资源采集站| 亚洲午夜高清在线观看| 人妻系列在线免费视频| 亚洲美女午夜激情视频在线观看| 综合激情网,激情五月| 国产av嗯嗯啊啊av| 国产精品美女免费视频观看| 婷婷六月天在线视频| 国产视频1区2区3区| 美国十次了亚洲天堂网国产| 国产视频成人一区二区| 日产国产欧美精品另类| 色噜噜噜噜色噜噜色合久一| 免费在线小视频你懂的| 成人av中文字幕在线看| 在线观看中文字幕精品av| 河北全程露脸对白自拍| 99精品久久99久久久久一| 天天干天天操天天日天天日| 羞羞漫画无限免费观看秋蝉| 精品免费一区二区三区四区视频| 抽插小穴啊啊啊视频| 可在线免费观看av| 国产igao激情在线视频入口| 一区二区三区国产在线成人av| 亚洲国产精品自拍偷拍视频在线| 九九视频在线观看全部| 中文字幕人妻一区色偷偷久久| 亚成区一区二区人妻熟女| 久久久国产精品免费视频网| 午夜福利片无码10000| 精产国品一二三77777| av男人站在线观看| av一区二区三区四区五区在线| 老牛影视在线一区二区三区| 深夜福利免费观看在线看| 日韩男女视频网站在线观看| 18禁男女啪啪啪无遮挡| 黑鸡巴肏少妇逼视频| 最新激情中文字幕视频| 91九色国产在线视频| 黑人3p日本女优中出| 国产精品 亚洲欧美 自拍偷拍| 最新久久这里只有精品| 男女爱爱好爽视频免费看| 欧美视频免费观看777| 欧美精品熟妇免费在线| 国产亚洲综合5388| 国产在线观看av一区| 亚洲男人天堂最新网址大全| 成人十欧美亚洲综合在线| 都市激情校园春色 亚洲| 日本少妇三级交换做爰做| 欧美vs亚洲vs日韩| 大尺度av毛片在线网址| 亚洲最大先锋资源采集站| 天天干夜夜爽狠狠操| 国产漂亮白嫩美女在线图片| 5d蜜桃臀女无痕裸感| 狠狠操深爱婷婷综合一区| 丰满少妇_区二区三区| 亚洲熟女乱色一区二区三区视频| 蜜乳av一区二区三区免费观看| 精品久久久久久久久久久久久 | 成人黄色录像在线观看| 久久无码高清免费视频| 人妻超清中文字幕在线乱码| 青青在线视频看看| 国产不卡免费在线观看| 最新久久这里只有精品| 青青草成人免费自拍视频| 九一精品人妻一区二区三区| 一区二区三区国产精华液区别大吗 | 亚洲免费午夜污福利| 天天操天天干天天谢| 人妻色综合aaaaaa网| 丝袜美女诱惑佐佐三上| 日本少妇熟女乱码一区二区| 国产夫妻视频在线观看免费| 日韩国产欧美久久一区| 少妇被中出一区二区| 一区二区三区四区视频精品免费| 中文字幕日本一二三区| 亚洲av中文无码网站| 四虎精品久久免费最新| 在线有码人妻自拍视频| 亚洲熟女乱色一区二区三区视频| 河北全程露脸对白自拍| av无限看熟女人妻另类av| 99久久国产精品免费热| 中文字幕熟女人妻一区| 女人高潮潮呻吟喷水网站| 午夜情色一区二区三区| 夏目彩春av在线看| 自拍偷自拍亚洲精品10p| 欧美大鸡吧男操女啊啊啊视频 | 午夜在线成人免费电影 | 在线看的免费网站黄| 国际日韩日韩日韩日韩日韩 | 最新中文字幕久久久久| 手机视频在线观看一区| 日本福利片在线播放| 午夜精品一区二区三区不卡顿| 啊不行啊操逼好爽大鸡吧视频| 自拍偷拍 国产激情| 日本不卡视频一二三区| 国产黑色丝袜 在线日韩欧美| 午夜精品久久秘?18免费观看| 黑人侵犯人妻森泽佳奈| 一区二区三区 国产日韩欧美| 十八禁黄色免费污污污亚洲| 天天曰天天摸天天爽| 黑人大巨屌操美女逼| av天堂hezyo| 在线免费观看a视频免费| 夜夜爽夜夜操夜夜爱| 精久久久久久久久久久久 | 一区二区三区四区久久久久韩日| 蜜臀一区二区日韩美女少妇视频 | 正在播放麻豆精品一区二区| 美女把逼扒开让男人桶| 中文字幕亚洲无线乱码| avtt中文字幕手机版| 亚洲一区二区三区四区入口| 亚洲国内精品久久久久久久| 国产精品免费看一区二区三区| 亚洲午夜国产末满十八岁勿进网站| 公侵犯人妻中文字幕巨| 日本电影一级人妻在线播放四区| 午夜精品视频免费观看| 日本高清久久人人爽| 成年人免费黄色av| 熟妇人妻av无码中文字幕| 日本人妻熟妇丰满成熟HD系列 | 美女av色播在线播放| 国产高清在线观看av| 91精品国产人妻麻豆| 日本人妻熟妇丰满成熟HD系列| 川上优所有中文字幕在线| 国产精品性感美女视频| 国产 亚洲 欧美 自拍| 精品精品精品精品精品污污污污 | 夜色福利视频免费观看| 日本有码精品一区二区三区| 中出小骚货在线观看| 欧美老熟妇xxoo老妇| 91精品国产成人久久久久久| 青青青青午夜手机国产视频| 中文乱码字幕人妻熟女人妻| 国产激情在线观看一区二区三区| 911美女片黄在线观看| 亚洲精品国品乱码久久久久| 91九色国产在线视频| 在线看日韩av不卡| 瑟瑟干视频在线观看| 日本不卡视频一二三区| 蜜乳av一区二区三区免费观看| 亚洲国产精品自产拍在线观看| 国产激情在线观看一区二区三区| 在线播放 日韩 av| 欧美日韩精品aaa| 欧美国产精品久久久免费| 免费看一级高潮喷水片| 亚洲国产电影的一区| 在线 激情 亚洲 视频| 99久久精品视频16| 日韩欧美中文字幕老司机三分钟| 可以免费观看日韩av| 日韩免费黄色片在线观看| 天天看片天天摸天天操| 精品精品精品精品精品污污污污| 最新免费在线观看污视频| 美女把逼扒开让男人桶| 亚洲av中文无码网站| 欧美一区二区三区爽爽| 美女把逼扒开让男人桶| 色就色综合偷拍区欧美在线| 亚洲一区二区三区无码在线| 18禁网站在线点击观看| 欧美亚洲另类精品第一页| 色欲AV蜜桃一区二区三| 不卡高清一区二区三区| 亚洲国产精品自拍偷拍视频在线| 日本高清在线观看不卡视频| 青娱乐免费最新视频| 大成色亚洲一二三区| 狠狠操深爱婷婷综合一区| 欧美黑人性猛交小矮人| 新香蕉视频香蕉视频2| av大尺度一区二区三区| 北野中文字幕一区二区| 天天夜夜久久精品综合| 亚洲理论在线a中文字幕97| 蜜桃臀av在线一区二区| 亚洲综合一区二区三区四区| 开心五月综合激情婷婷| 免费看一级高潮喷水片| tobu8日本高清| 久久人妻人人草人人爽| 天天操天天舔天天爽| 日韩精品视频一区二区三区在线| 亚洲国产日韩a在线欧美| 免费啪啪啪网站在线观看| 亚洲精品一区二区gif| 人妻系列级片在线观看视频| 成人av中文字幕在线看| 99热这里只有精品免费播放| 国产精品剧情av在线播放| 港台美女明星av天堂| 亚洲男人的天堂最新网址| 日本黄色一级电影网址| 极品风骚人妻3p视频| 欧美熟女xx00视频| 深夜福利免费观看在线看| 91porny九色视频偷拍| 超级黄肉动漫在线观看| 蜜乳av一区二区三区免费观看| 360偷拍蜜桃臀69式| 视频自拍偷拍视频自拍| 人人妻人人澡人人爽97| 伊人精品久久一区二区| 91亚洲最新蜜桃在线| 国产肥胖熟女又色又爽免费视频| 亚洲欧美国产一本综合首页| 欧美精品乱码99久久蜜桃免费| 天天日天天玩天天摸| 亚洲色大WWW永久网站| 99久9在线视频播放| 18禁网站在线点击观看| 欧美操大黑鸡巴视频在线观看| 欧美肥妇久久久久久| 婷婷综合缴情亚洲五月伊人| 不用付费特黄特色亚洲特级黄色片 | 成人午夜高清福利视频| 大香蕉在线欧美在线视频| 豆豆专区操逼性视频在线| 黑川堇人妻88av| 黑川堇人妻88av| 国产主播诱惑毛片av| 国产精品久久久久久成人久| 人妻系列在线免费视频| 国内自拍第一区二区三区| 夜夜人人干人人爱人人操| tushy一区二区三区视频| av天堂hezyo| 国产视频1区2区3区| 黄色av网址在线播放| 九九六视频,这里只有精品| 美女把腿张开给男的捅| 97香蕉久久国产超碰| 亚洲一区视频中文字幕在线播放| 成人18禁高潮片免费日本| 可在线免费观看av| 天天操天天干天天舔天天| 欧美黄色性视频网站| 丰满人妻熟女aⅴ一区| 松本菜奈实最新av在线 | 夜夜躁婷婷av蜜桃妖| 网站在线观看蜜臀91| 亚洲精品色图1234| 九九视频在线观看全部| 亚洲一区亚洲二区成人福利| 人妻少妇的va视频| 亚洲gay视频在线观看| 亚洲一区二区在线激情| 欧美日韩精品aaa| av网页免费在线观看| 中文字幕熟女乱一区二区| 亚洲熟女乱色一区二区三区视频| 4438x亚洲最大的成人| 午夜精品秘一区二区三区| 欧美日韩久久丝袜在线| 97精品国产91久久久| 北野中文字幕一区二区| 在线观看黄页网站视频网站| 漂亮人妻口爆久久精品| 婷婷综合缴情亚洲五月伊人| 91 精品视频在线看| 亚洲乱码国产乱码精品精视频| 久久一级片三上悠亚| jandara在线观看| 老鸭窝在线毛片观看免费播放| 美利坚合众国av天堂| 自拍偷拍亚洲综合第一页| 精品日本少妇久久久| 极品内射老女人操逼视频| 天天操天天舔天天射天天日天天干| 国产极品气质外围av| 精品高潮呻吟久久av| 亚洲国产精品青青草| 亚洲综合首页综合在线观看 | 日本老女人日比视频| 男人资源站中文字幕| 亚洲熟女乱色一区二区三区视频| 国产精品黄色片大全| 亚洲gay视频在线观看| 亚洲成人欧洲成人在线| 55夜色66夜色亚洲精品| 成人av在线视频免费| 最新日韩中文字幕免费在线观看 | 色视频免费观看网址| 国产一区两区三区福利小视频| 日本高清 中文字幕| 熟妇精品午夜久久久久| 国产精品久久人人添| 天堂av国产av伦理av| 亚洲av毛片一区二区三区网| 加勒比东京热绿帽人妻多人操| 日日夜夜免费视频精品| 亚洲一区二区中文字幕久久| 网站在线观看蜜臀91| 亚洲日本欧美韩国另类综合| 999久久久人妻精品一区| 亚洲少妇视频在线观看| 杜达雄啪啪毛片视频| 日本福利视频网站导航| av在线男人的天堂亚洲| 最近最新欧美日韩精品| 黄色av网址在线播放| 中文字幕熟女人妻一区| 999精品视频免费在线观看| 黑吊操欧美极品美女| 日本一区二区三区的资源| 国产精品福利久久久久| 欧美猛少妇色ⅹⅹⅹⅹⅹ猛叫| 亚洲午夜精品一级毛片app| 超级黄肉动漫在线观看| 超碰在线免费观看视频97| 在线观看中文字幕视频成人| 国产激情一区二区视频| 中文字字幕在线精品乱码| 亚洲另类欧美综合久久| 中文字幕熟女人妻丝袜丝在线| 92麻豆一区二区三区| 日本有码精品一区二区三区| 亚洲一区二区在线激情| yellow在线亚洲精品一区| 视频在线 一区二区| 天天在线播放日韩av| 中文字幕 首页 人妻| 综合久久伊人久久88| 可在线免费观看av| 欧美在线观看视频欧美| 久久无码高清免费视频| 国内精品一区二区2021在线| 香港日本台湾经典三级| 五十岁熟妇高潮喷水| 欧美成人性生活视频播放| a级片特黄免费看| 91精品91久久久久| 男生和女生羞羞91在线看| 豆豆专区操逼性视频在线| 97精品视频,全部免费| 国产自拍偷拍在线精品| 自拍偷拍亚洲综合第一页| 99热在线只有的精品| 午夜久久久久欠久久久久| 五月的婷婷综合视频| 中文字幕在线观看亚洲情色| 亚洲中文字幕无线乱码人妻精品 | 91人妻人人爽色啊啊啊| jizzjizz国产精品传媒| 亚洲av毛片一区二区三区网| 熟女人妻精品视频一区| 荣立三等功退休有什么待遇| 欧美成人短视频在线播放| 一区二区三区观看在线| 亚洲成人中文无码在线| 最近中文字幕免费视频一| 中文字幕 中文字幕 亚洲| 久久久久夜色国产精品电影| 日韩最近中文在线观看| 欧美亚洲精品色图网站| 人妻女侠被擒受辱记| 4438x亚洲最大的成人| 最近在线中文字幕免费| 九热精品视频在线观看| 天天摸天天舔天天操天天日| 午夜久久人妻一级内射av网址| 青青操91美女国产| 欧美黑人性猛交小矮人| 亚洲欧美不卡专业视频| 午夜精品秘一区二区三区| 熟妇高潮久久久久久久| 中文字幕观看中文字幕免费 | 欧美情色av在线观看| 啪啪啪网站免费在线看| 午夜92福利1000| 中字幕人妻熟女人妻a62v网| 亚洲精品9999蜜桃| 日本老女人日比视频| 一二区二区不卡视频| 国产清纯一区二区在线观看| 中文字幕观看中文字幕免费 | 人妻免费视频黄片在线视频| 亚洲中文字幕最新地址| 漂亮人妻口爆久久精品| 十八禁黄色免费污污污亚洲| 在线观看中文字幕视频成人| 天天干天天日天天弄| 欧美熟女xx00视频| 黑人大吊大战亚洲女人。| 熟妇高潮久久久久久久| 亚洲国产日韩a在线欧美| 欧美不卡一二三区精品| 国产男女无套?免费网站下载 | 99re这里是国产精品首页| 亚洲综合首页综合在线观看| 亚洲黄色成人一级片| 夜色福利视频免费观看| 中文字字幕在线精品乱码| 黄色av日韩在线观看| 亚洲精品激情视频在线观看| 日韩A级毛片免费视频| 自拍偷拍 亚洲性图 欧美另类| 欧美精品一区二区三区观看| 快色视频在线观看免费| 60路70路日本熟妇| 91精品一区一区三区| 18禁男女啪啪啪无遮挡| 九九热在线精品播放| 天堂网免费在线电影| 亚洲美女露隐私av一区二区精品| 丰满放荡熟妇在线播放| 91进入蜜桃臀在线播放| 丰满少妇人妻一区二区三区蜜桃 | 亚洲三级综合在线观看| 外国美女舔男人坤坤| 精品人妻在线激情视频| 大乳人妻一区二区三区| 91精品久久久久久久99蜜月 | 亚洲午夜精品一级毛片app| 国产av啊啊啊啊啊啊啊| 亚洲乱码av一区二区蜜桃av | 国产中文亚洲熟女日韩| 夜夜爽夜夜操夜夜爱| 91亚洲精品久久蜜桃| 免费看日韩黄视频在线观看| 亚洲精品激情视频在线观看| 天天在线播放日韩av| 免费观看在线中文字幕视频| 亚洲国产精品一区51动漫| 久久久久国产精品二区| 成人18禁高潮片免费日本| 美女福利视频一区二区三区四区| 性感美女人妻久久久| av里面的动作是真进去吗| 乱子伦国产一区二区三区| 少妇被中出一区二区| 97精品久久久久久无码人妻| 国产午夜在线播放视频| 亚洲一区二区三区国产精品电影| 插鸡视频免费网站在线播放| 美国十次了亚洲天堂网国产| 老熟女 露脸 嗷嗷叫| 老牛影视在线一区二区三区| 亚洲经典av中文字幕| 久久99国产中文丝袜| 69国产在线视频网站| 91大神福利视频网| 亚洲中文字幕无线乱码人妻精品| 亚洲美女a级黄色在线播放| 中字幕人妻熟女人妻a62v网| 每日更新日韩欧美在线| 免费在线观看亚洲福利| 午夜3p福利视频合集| 另类欧美激情校园春色| a级片特黄免费看| 99在线视频精品观看高| 亚洲熟妇在线视频观看| 日本高清有码在线视频| 中文字幕福利视频第四页| 久草视频在线看免费| 亚洲精品国品乱码久久久久| 中文字幕在线免费观看人妻| 秋霞成人午夜鲁丝一区二区三区| 一区二区三区不卡免费视频网站| 黑鸡巴肏少妇逼视频| 天天插天天干天天狠| 婷婷色综合五月天视频| 中文字幕 首页 人妻| 男人av一区二区三区| 欧美一级特黄大片做受99| 日本福利视频网站导航| 成年人免费福利在线| 亚洲午夜精品一级毛片app| 午夜精品小视频在线播放| 久久久久高潮白浆久久| 国产亚洲精品啪啪视频| 天天操天天干天天谢| 两个奶被揉得又硬又翘怎么回事| 最近日韩免费在线观看| 日韩国产欧美久久一区| 亚洲avav天堂av在线网毛片| 高潮喷水一区二区三区| 天天想要天天操天天干| 四季av人妻一区二区三区| 裸露视频免费在线观看| 免费在线观看亚洲福利| 亚洲无码专区中文字幕专区| 不卡一二三区别视频| 欧美区一区二区三视频| 快使劲弄我视频在线播放| 午夜夫妻性生活视频| 成人十欧美亚洲综合在线| 熟妇高潮久久久久久久 | 日本一道中文字幕99| 欧美成人区一区二区三| 亚洲一区二区在线激情| 日韩av熟妇在线观看| 久久久精品人妻无码专区不卡| 日本欧美国产在线一区| 中文在线字幕免费观看日韩视频| 亚洲熟女乱色一区二区三区视频| 免费24小时人妻视频| 不卡高清一区二区三区| 69精品互换人妻4p| 青娱乐免费最新视频| 亚洲无人区乱码中文字幕一区| 国语精品视频自产自拍| 蜜桃臀av在线一区二区| 人妻系列中文字幕大乳丰满人妻| 最新激情中文字幕视频| 一区二区三区 国产日韩欧美| 91亚洲精品久久蜜桃| 中文字幕人妻一区色偷偷久久| 黑人3p日本女优中出| 伊人网在线观看 视频一区| 大香蕉尹人在线最新| 国产熟女五十路一区二区三区| 亚洲国产中文字幕在线看| 蜜乳av中文字幕一区二区| 18岁禁一二三区免费体验| 亚洲中文字幕无线乱码人妻精品| 中字幕人妻熟女人妻a62v网| 亚洲欧美国产一本综合首页| 999久久久人妻精品一区| 久久99精品热在线观看| 性感人妻 中文字幕| 一区二区三区资源视频| 中文字幕在线观看av观看| 日韩一级视频一区二区三区| 红桃视频国产av在线| 日韩欧美一区二区三区免费看| 一区二区九日韩美女| 91在线九色porny| 18禁网站在线点击观看| 亚洲精品激情视频在线观看| 亚洲色大WWW永久网站| 美女福利视频一区二区三区四区| 亚洲国产精品一区51动漫| 成人资源中文在线观看| 亚洲人成大片在线观看| 九九热视频1这里只有精品| 午夜国产免费视频亚洲| 亚洲另类欧美综合久久| 新香蕉视频香蕉视频2| 亚洲一区二区偷拍女厕所| 91精品综合久久久久久五月天| 国内自拍第一区二区三区| 骚穴被阴茎插免费视频| 日本小视频一区二区| 在线观看中文字幕少妇av| 国产肥胖熟女又色又爽免费视频| 熟女国内精品一区二区三区| 亚洲 自拍 激情 另类| 亚洲国产精品一区51动漫| 69视频在线精品国自产拍| 天堂av国产av伦理av| 最近最新最好看的中文字幕| 99久久国产精品免费热| 99在线视频精品观看高| 欧美日韩久久丝袜在线| 国产精品性感美女视频| 亚洲成a人77777| 无人区一码二码三码区别在哪| 国产精品亚洲精品亚洲| 少妇精品视频一区二区免费看| 夜夜爽夜夜操夜夜爱| 国产午夜在线播放视频| 北野中文字幕一区二区| 高清国产美女a一级毛片| 精品国产久久久久午夜精品av| 人妻系列在线免费视频| 亚洲国产精品自拍偷拍视频在线| 国产精品剧情在线亚洲| 四虎国产精品国产精品国产精品| 性感美女极品18禁网站在线| 夜夜人人干人人爱人人操| 欧美成人久久久桃色aa| av 一区二区三区 熟女| 男人用大鸡巴狂操女人肉穴| 日本不卡视频一二三区| 久久久亚洲熟女一区二区| aaaa级少妇高潮在线观看 | 亚洲美女露隐私av一区二区精品 | 老司机免费视频福利0| 国产美女视频带a∨黄色片| av天堂hezyo| 亚洲成人动漫av在线| yy4080黄色片| 韩日一级人添人人澡人人妻精品| 一区二区三区不卡免费视频网站 | 港台美女明星av天堂| 午夜亚洲国产精品中字| 亚洲va999天堂va| 熟女人妻aⅴ一区二区三| 免费中文三级在线观看| 美利坚合众国av天堂| 日本韩国欧美在线视频| 可以免费观看日韩av| 久久久久久久久久久久久国产| av日韩视频在线观看| 亚洲综合熟女乱中文| 久久精品国产亚洲av热软件| 久草久热这里只有精品| 国语精品视频自产自拍| 91美女在线观看视频| 二十四小时日本高清在线观看 | 日本免费人爱做视频在线观看不卡 | 亚洲国产电影的一区| 亚洲精品综合欧美精品综合| 美女黄色啊啊啊啊视频| 黑人3p日本女优中出| 裸露视频免费在线观看| av 一区二区三区 熟女| 999国产精品视频免费看| 快使劲弄我视频在线播放| 最新国产精品久久精品app| 黑人和日本人av一区二区| 亚洲欧美精品日韩偷拍| 亚洲欧美激情国产综合久久久| 精品欧美黑人一区二区三区| 亚洲字幕一区二区夜色av| 西野翔人妻中文字幕中字在| 二十四小时日本高清在线观看| 99热99这里免费的精品| 日本有码精品一区二区三区| 日韩av熟妇在线观看| 99久久精品视频16| 欧美日韩综合精品无人区| av在线观看视频免费| 欧美成人屋影院在线视频观看| av成人三级高清日韩| 男人电影天堂在线观看| 黄在线看片免费人成视频| 夜夜操天天干夜夜操| 一区二区在线观看视频观看| 黄在线看片免费人成视频| 亚洲理论在线a中文字幕97 | 日韩人妻一区二区三区在线观看| 国产主播诱惑毛片av| 伊人久久综合国产精品| 男人电影天堂在线观看| 大香蕉在线欧美在线视频| 亚洲综合一区二区三区四区| 77亚洲视频在线观看| 特级aaaaa黄色片| 猫咪亚洲中文在线中文字幕| 日本韩国欧美在线视频| 鸡巴插进美女的嫩小穴视频| 国产肥胖熟女又色又爽免费视频 | av无限看熟女人妻另类av| 国产毛片特级Av片| 在线观看中文字幕少妇av| 黄色av 在线观看| 精品人妻在线激情视频| 亚洲AV无码一二三四区在线播放| 亚洲 自拍 激情 另类| 欧美色区国产日韩亚洲区| 熟女人妻少妇一区二区| 两个人在一起靠逼啊啊啊| 国产乱码有码一区二区三区| 蜜臀一区二区日韩美女少妇视频| 国产激情免费在线视频| 精品不卡一区二区三区| 久久内射天天玩天天懂色| 天天天天天天天天干夜夜| 美利坚合众国av天堂| 国产青青青青草免费在线视频| 日韩加勒比精品在线看| 夜夜人人干人人爱人人操| 亚洲免费在线不卡视频| 精品国模一区二区三区欧美| av一区二区三区蜜桃| 亚洲经典av中文字幕| 日本国产亚洲欧美色综合| 汤姆提醒30秒中转进站口| 亚洲色大WWW永久网站| 天天操天天舔天天射天天日天天干| 天天做天天日天天搞| 人妻系列中文字幕大乳丰满人妻| 青青青免费手机视频在线观看| 亚洲国产日韩a在线欧美| 久久久久久a女人处女| 美利坚合众国av天堂| 夜夜爽夜夜操夜夜爱| 色丁香久久激情综合网| 亚洲美女午夜激情视频在线观看 | 69久久夜色精品国产69乱电影| 白白色在线免费视频发布视频| 中文人妻av一区二区三区| 亚洲欧美韩国日本一区二区| 色就色综合偷拍区欧美在线| 91精品久久久久久久99蜜月| 国产 少妇 一区二区| 麻豆出品视频在线观看| 久久久国产精品免费视频网| ass亚洲熟女ass| 美女欧美视频在线观看免费| 欧美vr专区日韩vr专区| 欧美亚洲精品色图网站| 91精品国产成人久久久久久| 91久久久久久最新网站| 开心五月综合激情婷婷| 欧美啪啪一区二区三区| 麻豆午夜激情在线观看| 亚洲AV无码久久精品国产一区老| 欧美日韩亚洲tv不卡久久| 国产精品美女免费视频观看| 夫亡人妻被强干中文字幕| 人人妻人人爽人人摸| 日本黄色一级电影网址| 美女激情久久久久久久| 免费中文三级在线观看| 国产资源在线观看二区| 成人午夜av电影网| 新香蕉视频香蕉视频2| 中文字幕亚洲乱码精品无限| 亚洲精品色图1234| 精品一区二区三区喷水内射高潮 | av日韩视频在线观看| 亚洲欧美一级特黄大片 | 人妻色综合aaaaaa网| 国产经典精品欧美日韩| 青青青在线视频观看97| 亚洲黄色成人一级片| 岳母的诱惑电影在线观看| 亚洲午夜熟女在线观看| 69视频在线精品国自产拍| 岳母的诱惑电影在线观看| 18岁禁一二三区免费体验| av日韩视频在线观看| av在线播放观看h| 丰满放荡熟妇在线播放| 黄色av日韩在线观看| 中文字幕日韩首页欧美在线激情 | 18在线观看免费观看| 中文字幕观看中文字幕免费| 激情久久在线免费观看视频| 五月婷婷伊人久久中文字幕| 亚洲全国精品女人久久久| 欧美视频亚洲视频在线| www国产亚洲精品久久久| 裸日本资源在线午夜| 亚洲图片另类综合小说| 日韩国产欧美久久一区| 大香蕉伊人97在线| 国产白丝一区二区三区av| 青青操91美女国产| 天天看天天爱天天日| 最近中文字幕免费视频一| av成人三级高清日韩| 欧美日本国产一区二区| 亚州av嫩草av极品在线观看| tushy一区二区三区视频| 精品一区二区三区免费毛片W| 四虎国产精品国产精品国产精品| 天天弄天天草天天日天天| 人妻免费视频黄片在线视频| 亚洲综合色一区二区三区| 99精品久久一区二区| 国产视频成人一区二区| xxoo福利视频导航| 中文字幕日韩人妻在线三区| 中文字幕综合网91| 人妻超清中文字幕在线乱码| 亚洲精品激情视频在线观看| 黄色片免费国产精品| 国产在线观看av一区| 亚洲欧美另类丝袜另类自拍| 国产精品剧情在线亚洲| 日本五六十路熟女视频| 亚洲成年人精品国产| 干逼又爽又黄又免费的视频| 在线有码人妻自拍视频| 在线观看网站伊人网| 丰满放荡熟妇在线播放| 日本高清久久人人爽| —区二区三区女厕偷拍| 啊~插得好快别揉我胸了视频| 2021国产在线视频| 国内精品一区二区2021在线| 久久久人妻免费视频| 午夜久久人妻一级内射av网址 | 人人妻人人狠人人爽| 内地精品毛片在线观看| 国产做A爱免费视频在线观看| 国产视频成人自拍蝌蚪视频 | 女女抠逼白虎白丝袜| 亚洲av中文无码网站| 天天爱天天日天天爽| 亚洲另类激情视频在线看| 亚洲欧美成人激情在线| 欧美熟女xx00视频| 亚洲一区在线视频观看地址| 夜夜爽夜夜操夜夜爱| 亚洲制服丝袜网站中文字幕| 网友自拍第一页99热| 久久久久性感美女偷拍视频| 中文字幕福利视频在线一区| 天天摸天天干夜夜操| 2018中文字字幕人妻| 午夜美女福利视频在线| 中文字幕免费啪啪啪| 亚洲蜜桃久久久久久| 一区二区三区五区六区| 日本福利视频网站导航| 人妻激情偷乱一区二区三区av| 天天夜夜久久精品综合| 亚洲|久久久久久一二三区丝袜| 日本不卡 中文字幕| 成人精品动漫一区二区| av在线免费在线观看| 亚洲欧美另类丝袜另类自拍| 天天爽天天操天天插| yellow在线亚洲精品一区| 亚洲美女色www色| 自拍偷自拍亚洲精品10p| 在线成人教育平台排名| 中文字幕丰满子伦无码专区| 国产亚洲精品啪啪视频| 大成色亚洲一二三区| 国产资源网站在线播放| 韩国一级片最火爆中文字幕| 中文字幕熟女人妻一区| 色视频在线播放免费观看| 国产,亚洲,欧美综合| 老熟妇一区二区三区v∧88| 午夜久久人妻一级内射av网址| 两个人在一起靠逼啊啊啊| 伊人免费观看视频一| 成人黄色录像在线观看| 亚洲午夜熟女在线观看| 人人妻人人爽人人摸| 国产激情一区二区视频| 中文字幕熟女乱一区二区| 激情久久在线免费观看视频| 人妻超清中文字幕在线乱码| 一区二区三区五区六区| 成人黄色录像在线观看| 亚洲av综合av一去二区三区| 中文字幕熟女乱一区二区| 丰满人妻熟女aⅴ一区| 黑人巨大精品一区二区在线| 91污污在线观看视频| 中国精品人妻一区二区| 久久久精品人妻无码专区不卡| 午夜野花视频在线观看| 日韩A级毛片免费视频| 韩国资源视频一区二区三区 | 午夜精品久久久久久久精品乱码| 国产美女视频带a∨黄色片| 女同性恋av在线播放| 免费在线观看黄色小网站| 快色视频在线观看免费| 一级做性色a爱片久久片| 亚洲人成小说网站色| 人妻超清中文字幕在线乱码| 在线免费视频999| 精品精品精品精品精品污污污污| 91激情四射婷婷综合| 国产黑色丝袜 在线日韩欧美| 久久久久夜色国产精品电影| 亚洲成人中文无码在线| 久久中文字幕av一区二区| 91精品一区一区三区| 婷婷六月天在线视频| 中文字幕人妻一区色偷偷久久| 51vv精品视频在线观看| 成人av中文字幕在线看| 亚洲自拍偷拍av在线| 女人扒开逼让男人操| 人妻被强av系列一区二区| 亚洲黄色成人一级片| 国产igao激情在线视频入口| 五十岁熟妇高潮喷水| 日本少妇人妻凌辱在线| 青青在线视频看看| 国产成人av在线你懂得| 欧美日韩亚洲tv不卡久久| 911美女片黄在线观看| 在线看日韩av不卡| 国产av剧变态维修工虐杀美女| 4438全国成人免费视频| 久久视频 在线播放| 日本高清激情乱一区二区三区| 午夜在线观看一级毛| 一区二区三区观看在线| 亚洲免费午夜污福利| 老司机伊人99久久精品| 天天插天天操天天射天天干| 91精品在线视频免费视频| 日韩女同与成人用品电影免费看| 日本欧美国产在线一区| 丰满人妻被猛烈进入中文字幕| 奇米网首页神马久久| 极品少妇高潮喷水日出白浆| 亚洲第一页欧美第一页| iga肾三级算严重吗| 欧美男男在线观看视频网站| 久久人妻诱惑我视频| 日本五六十路熟女视频| 最新日韩中文字幕啪啪啪| 精品国模一区二区三区欧美| 天天日夜夜操人人爽| 中文字幕熟女人妻一区| 欧美肥妇久久久久久| 亚洲一区二区偷拍女厕所| 亚洲国产精品自拍偷拍视频在线| 91超碰九色porny| av日韩视频在线观看| 日韩三级黄色大片在线观看| 最新久久这里只有精品| 自拍偷自拍亚洲精品10p| 女女抠逼白虎白丝袜| 欧美一级特黄大片在线| 久久免费视频ww一区| 久久国产精品久精国产爱| 5d蜜桃臀女无痕裸感| 亚洲成人偷拍自拍在线| 久久av色噜噜ai换脸| 2021国产剧情麻豆| 天天干天天色综合久久| 午夜在线成人免费电影| 亚洲国产精品 久久久| 九九六视频,这里只有精品| 欧美日韩在线观看免费播放| 午夜美女福利视频在线| 欧美成人短视频在线播放| 欧美视频免费观看777| 夜夜骚av一二三区| 日本电影一级人妻在线播放四区| 人妻少妇视频系列视频在线 | 欧洲成熟女人色惰片| 欧美性感美女热舞视频| 国产av精品一区二区三区久久| 天堂av国产av伦理av| 伊人免费观看视频一| 极品内射老女人操逼视频| 99re这里是国产精品首页 | 亚洲成人自拍av在线| 在线观看中文字幕少妇av| 亚洲熟女乱一区二区精品成人| 欧美日韩综合精品无人区| 激情久久在线免费观看视频| 亚洲|久久久久久一二三区丝袜| 自拍偷拍亚洲综合第一页| 色欲AV蜜桃一区二区三| 韩国资源视频一区二区三区| 亚洲av中文免费在线| jizzjizz国产精品传媒| 日韩国产欧美久久一区| 国产av嗯嗯啊啊av| 一区二区三区婷婷中文字幕| 亚洲成人av在线一区二区| 亚洲成人欧洲成人在线| 日本成人福利电影网| 三级欧美日韩一区二区三区| 九九九九九久久久国产| 好看的日本中文字幕在线观看二区| 亚洲av三级电影在线观看| 人人妻人人爽人人摸| 人妻被强av系列一区二区| 一区二区三区高清视频3| 伊人网在线免费观看| 老司机免费视频福利0| 日韩人妻中文字幕区| 中文字幕观看中文字幕免费| 人妻熟女 亚洲 一页二页| 高潮喷水一区二区三区| 亚洲午夜高清在线观看| 网站在线观看蜜臀91| 99国产精品国产精品毛片19| 日韩A级毛片免费视频| 国产亚洲综合5388| 国产福利一区二区三区在线观看| av网页免费在线观看| 美女精品久久久久久久久| 一区二区三区午夜福利在线| 青青操久久综合激情| 亚洲精品中文字幕手机在线免费看| 国产熟妇色xxⅹ交白浆视频| 国产肥胖熟女又色又爽免费视频 | 久久久久久久久久久久久国产| 人妻激情综合久久久久蜜桃| 欧美丝袜亚洲国产日韩| 68视频在线免费观看| 污视频在线观看地址| 久久内射天天玩天天懂色| 豆豆专区操逼性视频在线| 无人区一码二码三码区别在哪| 大尺度av毛片在线网址| 亚洲一区二区精品在线播放| 熟女俱乐部jukujoclub| 麻豆出品视频在线观看| 夜夜人人干人人爱人人操| av人摸人人人澡人人超碰小说| 男人电影天堂在线观看| 日韩欧美国产一区二区在线观看| 蜜臀一区二区日韩美女少妇视频| 蜜臀一区二区日韩美女少妇视频| 在线 制服 中文字幕 日韩| 亚洲精品9999蜜桃| 五十岁熟妇高潮喷水| 豆豆专区操逼性视频在线| 高清欧美色欧美综合网站| 正在播放麻豆精品一区二区| 18福利视频在线观看| av在线免费在线观看| 中文字幕av特黄毛片| 人妻色综合aaaaaa网| 操死你美女在线视频| 蜜桃臀少妇白色紧身裤细高跟| 国产成人在线观看视频播放| 国产精品亚洲精品亚洲| 综合激情网,激情五月| 手机看电影一区二区三区| 68福利精品在线视频| 一区二区三区国产精华液区别大吗| 日本电影一级人妻在线播放四区 | 亚洲制服丝袜在线看| 在线看的免费网站黄| 天天透天天舔天天操| 台湾18禁久久久久久久激情视频| 蜜桃tv一区二区三区| 成年男女免费视频网站无毒| 日本欧美国产在线一区| 一区二区三区四区 在线播放| 欧美国产精品久久久免费| 亚洲国产日韩a在线欧美| 自拍偷拍色图亚洲天堂| 国产午夜在线播放视频| 二十四小时日本高清在线观看| 国语对白性爱三级片免费看| 亚洲少妇视频在线观看| 男人的天堂av中文字幕| 1级黄色片在线观看| 精品一区二区三区免费毛片W| 亚洲黑人欧美二区三区| 一区二区三区内射美女| 漂亮人妻口爆久久精品| 免费在线小视频你懂的| 欧美性受黑人猛交裸体视频| 成人午夜高清福利视频| 在线观看网站伊人网| 天天透天天舔天天操| 久久精品四虎夜夜拍拍拍| 亚洲午夜高清在线观看| 精品国产人伦一区二区三区| 亚洲精品一区二区gif| 自拍偷拍 国产激情| 中文字幕人妻一区二区视频系列| 黑人大吊大战亚洲女人。| 亚洲天堂av最新在线| 豆豆专区操逼性视频在线| 免费中文字幕a级激情| 天堂网成人av电影| 久久久亚洲熟女一区二区| 伊人网国产在线播放| 国产av精品一区二区三区久久| 亚洲va999天堂va| 可以免费观看日韩av| 51vv精品视频在线观看| 最近最新最好看的中文字幕| 男人和女人的逼视频| 欧洲亚洲一区二区三区四区| 精品日本少妇久久久| 插鸡视频免费网站在线播放| 国产熟妇色xxⅹ交白浆视频| 一区二区在线观看视频网站 | 精品人妻在线激情视频| 欧美成人红桃视频在线观看| 国产午夜在线播放视频| 欧美巨大另类极品video| 美女黄色啊啊啊啊视频| 亚洲一区二区三区四区入口| 美女露阴道让男人捅| 亚洲国产综合久久精品| 青青青在线视频免费播放| 日本亚洲精品视频在线观看| 熟女人妻aⅴ一区二区三| 中文字幕在线免费观看人妻| 亚洲春色av中文字幕| av天堂a亚洲va天堂va里番| 97视频人人爱麻豆| 99久久99九九九99九| 午夜久久久久久av五月| 九九热在线精品播放| 白白色在线免费视频发布视频| 少妇熟女天堂网av| 午夜偷拍的视频久久久免费大全| 中文字幕人妻精品精品| 69精品人妻久久久久久久久久久| 色网站在线观看免费| 人人妻人人爽人人爽欧美一区| 亚洲少妇视频在线观看| 黑人巨大精品一区二区在线| 中文字幕一区二区三区久久久| 欧美一区二区三区视频看| 青青青在线视频观看97| 亚洲永远av在线播放| 嗯~嗯~啊啊啊~高潮了软件| 十八禁黄色免费污污污亚洲| 国产视频成人一区二区| 日本福利视频网站导航| 999精品视频免费在线观看| 日本亚洲午夜福利一区二区三区| 亚洲一区在线视频观看地址| 天天日 天天舔 天天射| 91精品一区一区三区| 亚洲午夜熟女在线观看| 欧美成人屋影院在线视频观看| 啪啪啪网站免费看视频| 亚洲人成大片在线观看| 日本香港韩国三级黄色| 欧美日本国产一区二区| 欧美精品999不卡| 久久久精品人妻无码专区不卡 | 加勒比不卡在线视频| 在线免费观看a视频免费| 丰满放荡熟妇在线播放| 五十岁熟女高潮喷水| 人妻在线中文视频视频| 日本成人福利电影网| 中文字幕熟女乱一区二区| 国产肥胖熟女又色又爽免费视频| 亚洲精品乱码久久久久app | 在线 激情 亚洲 视频| 日本清纯中文字幕版| 黑人侵犯人妻森泽佳奈| 日韩美精品成人一区二区三区四区| 美女张开腿给男人桶爽的软件| 在线中文字幕人妻av| 麻豆国产91制片厂| 亚洲人人爽人人澡起碰av| 97精品国产91久久久| 51vv精品视频在线观看| 欧美一区二区三区视频看| 女生抠逼自慰啊啊啊啊啊啊啊下载| 最新中文字幕久久久久| 亚洲制服丝袜资源网| 亚av一二三在线观看| 精品精品精品精品精品污污污污| 91国产精品乱码久久久久久| 午夜在线成人免费电影 | 日本熟妇乱妇熟色视频| 国产黄色主播网址大全在线播放| 性色蜜桃臀x88av天美传媒| 亚洲中文字幕在线视频观看二区| 久久中文字幕av一区二区| 五月婷婷伊人久久中文字幕| 日本香港韩国三级黄色| 亚洲欧洲无码一区2区无码| 夜夜操天天干夜夜操| 亚洲欧美另类丝袜另类自拍| 欧美激情视频第一页| 91系列视频在线播放| 天天日天天亲天天操| 人妻视频网站快射视频网站| 免费24小时人妻视频| 国产成人深夜福利短视频99| 狠狠操深爱婷婷综合一区| 亚洲av激情综合网| 99色在线观看免费观看| 狂操鸡巴小骚逼视频免费观看| 老牛影视在线一区二区三区| 欧美情色av在线观看| 美女张开腿给男人桶爽的软件| 一级毛片特级毛片免费的| 日韩成人在线电影首页| 黄片操操操操操操c| 亚洲午夜国产末满十八岁勿进网站 | 亚洲激情噜噜噜久久久| 欧美插插插插插插| 亚洲欧美精品海量播放| 国产农村乱子伦精精品视频| 99久久久久久久久久久久久| 熟女俱乐部jukujoclub| iga肾三级算严重吗| 天天日天天干天天日天天干天天 | 免费中文字幕a级激情| 天天色天天射天天日天天干| 免费观看在线中文字幕视频| 97人妻av人人澡人人爽| 亚洲成人 国产精品| 欧美精品一区二区三区观看| 91福利高清在线播放| 97人妻在线视频自拍| 97精品视频,全部免费|