水的氫鍵網(wǎng)絡(luò)動(dòng)力學(xué)與其太赫茲頻譜的關(guān)系*

段銅川 閆韶健 趙妍 孫庭鈺 李陽(yáng)梅 朱智?

1) (上海理工大學(xué)光電信息與計(jì)算機(jī)工程學(xué)院，教育部醫(yī)用光學(xué)儀器與設(shè)備重點(diǎn)實(shí)驗(yàn)室，上海 200093)

2) (國(guó)防科技創(chuàng)新研究院，太赫茲生物物理創(chuàng)新工作站，北京 100071)

水是萬(wàn)物生命之源，認(rèn)識(shí)水的太赫茲吸收譜是太赫茲技術(shù)在生物醫(yī)學(xué)上應(yīng)用的前提，太赫茲頻率的選擇對(duì)高效、低能耗地實(shí)現(xiàn)太赫茲的生物效應(yīng)至關(guān)重要.水的復(fù)雜氫鍵網(wǎng)絡(luò)使得其具有較寬的太赫茲吸收峰，因此有必要研究水的氫鍵網(wǎng)絡(luò)動(dòng)力學(xué)與其太赫茲吸收譜之間的關(guān)系，然而這方面的研究仍然非常缺乏.采用分子動(dòng)力學(xué)模擬方法，本文研究了不同水模型在常溫常壓下的太赫茲吸收譜，并且進(jìn)一步基于溫度研究了水的太赫茲吸收譜對(duì)氫鍵網(wǎng)絡(luò)強(qiáng)弱的依賴性，發(fā)現(xiàn)溫度的升高會(huì)使氫鍵網(wǎng)絡(luò)的太赫茲吸收譜發(fā)生紅移，這表明氫鍵網(wǎng)絡(luò)的太赫茲吸收譜的中心頻率與氫鍵相互作用的強(qiáng)弱具有強(qiáng)關(guān)聯(lián)，更進(jìn)一步的研究表明水中氫鍵網(wǎng)絡(luò)的氫鍵壽命與氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰的中心頻率之間存在線性關(guān)系.這一現(xiàn)象背后的物理能夠通過將氫鍵網(wǎng)絡(luò)中的氫鍵類比為彈簧借助彈簧振子模型加以描述.本文的發(fā)現(xiàn)將有利于理解水中復(fù)雜的氫鍵網(wǎng)絡(luò)動(dòng)力學(xué)，以及促進(jìn)太赫茲的生物效應(yīng)研究.

1 引言

水是生命之母[1].水在分子細(xì)胞生物學(xué)中表現(xiàn)出不同的結(jié)構(gòu)和功能，它作為溶劑幫助細(xì)胞內(nèi)的化學(xué)和信息傳遞過程，決定并參與生物分子的相互作用和運(yùn)動(dòng)[2-7].因此，研究水的動(dòng)力學(xué)性質(zhì)對(duì)于理解生化反應(yīng)以及大分子的生物功能等科學(xué)問題至關(guān)重要[8-12].隨著計(jì)算機(jī)的普及和計(jì)算能力的提升，使用計(jì)算機(jī)進(jìn)行分子動(dòng)力學(xué)模擬已經(jīng)成為與實(shí)驗(yàn)研究平行的一種研究方法.1971 年Rahman 和Stillinger[13]首次基于計(jì)算機(jī)模擬研究了具有分子團(tuán)簇行為的水的性質(zhì).為研究生物分子在溶液中的相互作用，研究人員提出了許多水模型.1981 年，Berendsen 等[14]基于簡(jiǎn)單性的要求提出了SPC 水模型，設(shè)計(jì)了一種新的液態(tài)水的有效對(duì)勢(shì)，采用由盡可能少的點(diǎn)電荷組成的模型，在保證足夠的計(jì)算精度的前提下，使得分子動(dòng)力學(xué)模擬研究體相水的動(dòng)力學(xué)性質(zhì)成為了可能.盡管該模型在大多數(shù)情況下的表現(xiàn)令人滿意，但仍有改進(jìn)的余地.1983 年，Jorgensen 等[15]提出了TIP3P 和TIP4P 水模型，TIP3P 能夠很好地描述生物分子在水溶液中的動(dòng)力學(xué)性質(zhì)，因此常用于生物系統(tǒng)[16]，TIP4P 是四點(diǎn)水模型，負(fù)電荷位于氫氧氫平分線上距離氧原子0.015 nm 的虛原子，而非氧原子上，這有效地改善了水分子周圍的靜電分布，它也常用于生物體系.1987 年，Berendsen 等[17]對(duì)SPC 模型進(jìn)行了重新參數(shù)化，以獲得正確的密度和能量，由此得到了SPC/E水模型.2004 年Horn 等[18]在TIP4P 的基礎(chǔ)上又提出了TIP4P-Ew 水模型，被廣泛使用于Ewald 求和方法.盡管如此，Guillot[19]對(duì)不同水模型進(jìn)行了深入的研究，發(fā)現(xiàn)特定水模型只能反映水在某一方面的性質(zhì)，迄今沒有一個(gè)令人完全滿意的水勢(shì)模型.

最近的研究表明，水的動(dòng)力學(xué)過程與太赫茲之間有著密切的關(guān)系.太赫茲(terahertz，THz) 波的波段能夠覆蓋生物大分子、有機(jī)體、半導(dǎo)體和等離子體等物質(zhì)的特征譜，利用該頻段可以加深和拓展人類對(duì)物理學(xué)、化學(xué)、天文學(xué)、信息學(xué)和生命科學(xué)中一些基本科學(xué)問題的認(rèn)識(shí)[20-24].在生物醫(yī)學(xué)上，生物大分子相互作用是重大生命現(xiàn)象與病變產(chǎn)生的關(guān)鍵動(dòng)因，然而太赫茲光子的能量覆蓋了生物大分子空間構(gòu)象的能級(jí)范圍 (如圖1(a) 所示)，該頻段包含了其他電磁波段無(wú)法探測(cè)到的直接代表生物大分子功能的空間構(gòu)象等重要信息[25].2018 年，劉國(guó)治院士[26]推斷生物神經(jīng)信號(hào)物理場(chǎng)應(yīng)為太赫茲到紅外區(qū)域的高頻電磁場(chǎng)，最可能頻率范圍應(yīng)在0.5—100 THz，并且將此波段稱為廣義太赫茲.這預(yù)示著廣義太赫茲波有望調(diào)控生物分子的結(jié)構(gòu)和功能.生物分子在水溶液中發(fā)揮活性和功能，然而水是極性液體，它在部分太赫茲頻率下有異常高的吸收損失，表現(xiàn)出熱效應(yīng)，給太赫茲技術(shù)在生物體系中的應(yīng)用蒙上了一層烏云[27].通常來(lái)說，水的太赫茲光譜有3 個(gè)主峰 (圖1(b))，第1 個(gè)峰位于5—30 THz 區(qū)間，對(duì)應(yīng)水中氫鍵網(wǎng)絡(luò)的振動(dòng)模式(圖1(c))，第2 個(gè)峰位于45—50 THz，對(duì)應(yīng)水分子內(nèi)部鍵角的彎曲振動(dòng) (圖1(c))，第3 個(gè)峰則位于90—105 THz，它對(duì)應(yīng)于水分子內(nèi)部鍵長(zhǎng)的拉伸振動(dòng) (圖1(c)).這意味著在廣義太赫茲波段存在4 個(gè)水的太赫茲弱吸收窗口，可用于非熱地調(diào)控生物大分子的結(jié)構(gòu)和功能.得益于此，近3 年，太赫茲在水溶液中非熱地調(diào)控細(xì)胞動(dòng)力學(xué)過程取得重要的進(jìn)展[28-35](圖1(d)).2019 年，我們提出1.39 或4.66 THz 的電磁刺激能夠非熱地引發(fā)一維水通道內(nèi)受限水的超級(jí)滲透，促使受限水由一維冰相向一維相干氣相轉(zhuǎn)變[28].2020 年，Wu 等[29]發(fā)現(xiàn)44.0 THz 的電磁刺激能夠共振地加速DNA 分子的解旋過程，并且有效地降低DNA 的熔解溫度.同年，Wang 等[30]觀察到太赫茲輻射能夠引起DNA堿基對(duì)中氫鍵的瞬態(tài)質(zhì)子轉(zhuǎn)移.Li 等[31]證實(shí)細(xì)胞內(nèi)ATP 水解能夠釋放頻率約為34 THz 的內(nèi)源光子，進(jìn)而調(diào)控生化反應(yīng).2021 年，Li 等[32]發(fā)現(xiàn)42.55 THz 的電磁刺激能夠非熱地加速鈣離子通道內(nèi)鈣離子的滲透.同年，Liu 等[33]演示了53.7 THz的電磁刺激能夠非熱、可逆地調(diào)控神經(jīng)信號(hào)和動(dòng)物的行為.以上研究表明，為了避免太赫茲電磁波被體相水強(qiáng)吸收，理論上預(yù)測(cè)既能夠調(diào)控生物大分子的功能，又對(duì)體相水產(chǎn)生有限熱效應(yīng)的太赫茲波段非常關(guān)鍵，因此，有必要深入地對(duì)體相水的太赫茲吸收譜進(jìn)行理論研究.

圖1 太赫茲與生物分子的密切關(guān)系以及太赫茲調(diào)控細(xì)胞動(dòng)力學(xué) (a) 生物大分子的轉(zhuǎn)/振動(dòng)的頻率在THz 頻段;(b) 水的太赫茲吸收譜，綠色區(qū)域是有望非熱地調(diào)控生物分子的廣義太赫茲頻率的4 個(gè)窗口;(c) 水的太赫茲吸收譜的振動(dòng)模式[36];(d) 太赫茲波非熱調(diào)控細(xì)胞動(dòng)力學(xué)，涉及發(fā)揮細(xì)胞生物功能的水通道蛋白、DNA、鉀離子通道、鈣離子通道Fig.1.Close relationship between terahertz and biomolecules and the regulation of cell dynamics by terahertz:(a) Frequency of rotation/vibration of biological macromolecules is in the THz frequency band;(b) terahertz absorption spectrum of water，the green region is the four frequency windows in which electromagnetic wave is expected to non-thermally regulate biomolecules;(c) vibration modes of water corresponding to its terahertz absorption spectrum[36];(d) terahertz waves non-thermally regulate the dynamics of a cell，involving aquaporins，DNA，potassium and calcium channels that perform biological functions of the cell.

本文利用分子動(dòng)力學(xué)模擬方法基于不同水模型和溫度，研究了水的太赫茲吸收譜以及與氫鍵網(wǎng)絡(luò)動(dòng)力學(xué)的關(guān)系，發(fā)現(xiàn)水形成的氫鍵網(wǎng)絡(luò)越強(qiáng)，其太赫茲吸收峰越傾向于藍(lán)移，進(jìn)一步發(fā)現(xiàn)了水的氫鍵壽命與其太赫茲吸收峰的頻率之間存在線性關(guān)系，并利用彈簧振子模型對(duì)這一現(xiàn)象的機(jī)理加以解釋，也就是將水的氫鍵等同于彈簧，如果水的氫鍵越強(qiáng)，其等效彈簧的勁度系數(shù)就越大，進(jìn)而使得氫鍵網(wǎng)絡(luò)具有較大的振動(dòng)頻率.

2 研究方法

2.1 模擬方法

首先，建立大小為4 nm × 4 nm × 4 nm 的立方體模擬盒子，然后在常溫常壓的情況下分別用水分子個(gè)數(shù)為2165 的柔性SPC，SPC/E，TIP3P 水分子模型和水分子個(gè)數(shù)為2211 的柔性的TIP4P，TIP4P-Ew 水分子模型將模擬盒子填充滿.水分子間的相互作用為

式中第一項(xiàng)為長(zhǎng)程靜電相互作用，第二項(xiàng)為短程Lennard-Jones作用.其中f1/(4πε0)，ε0是真空中的介電常數(shù);rij是原子i和j之間的距離;qi是原子i的電荷;εij，σij是原子i和原子j之間的Lennard-Jones 作用參數(shù).這些水模型的力場(chǎng)參數(shù)差異如表1 所列，其中，rH—O，θ分別表示O—H 的鍵長(zhǎng)、∠ HOH 的鍵角，qO(qvir)，qH為氧原子(虛原子)、氫原子所帶電荷量，ε和σ為L(zhǎng)ennard-Jones作用參數(shù).模擬使用的軟件是GROMACS 5.1.4，采用OPLS/AA 力場(chǎng)處理原子間的相互作用，其中不同原子間的Lennard-Jones 參數(shù)基于幾何平均法計(jì)算，使用跳蛙算法求解粒子的運(yùn)動(dòng)方程.模擬過程采用周期性邊界條件，模擬系綜為恒溫恒壓系綜，使用Nose-Hoover 熱浴控制系統(tǒng)溫度保持在設(shè)定的參考溫度，溫度弛豫時(shí)間為0.2 ps，使用Parrinello-Rahman 壓強(qiáng)耦合器控制系統(tǒng)壓強(qiáng)維持在1 bar (1 bar=105Pa).使用Particle-Mesh-Ewald 方法處理靜電相互作用，靜電和范德瓦耳斯相互作用的截?cái)喟霃骄鶠? nm.每次模擬時(shí)長(zhǎng)為100 ns，時(shí)間步長(zhǎng)為1 fs，其中最后5 ns 的軌跡數(shù)據(jù)用于分析水的光譜和動(dòng)力學(xué)性質(zhì).

表1 所用水模型的力場(chǎng)參數(shù)的比較Table 1. Comparison of the force field parameters of different water models employed.

2.2 分析方法

2.2.1 頻譜解析

在經(jīng)典近似的情況下，分子與頻率相關(guān)的光譜強(qiáng)度(I)可以通過對(duì)水的總電荷流的自相關(guān)函數(shù)進(jìn)行傅里葉變換獲得，它的表達(dá)式為

2.2.2 徑向分布函數(shù)

徑向分布函數(shù)表示了與某個(gè)原子距離為r的單位體積元內(nèi)出現(xiàn)的分子數(shù)密度，表達(dá)式為[39，40]

其中，N是系統(tǒng)中水分子的個(gè)數(shù)，ρN是水的平均密度，rij是第i號(hào)分子和第j號(hào)分子之間的距離.

2.2.3 自擴(kuò)散系數(shù)

根據(jù)Einstein 關(guān)系，自擴(kuò)散系數(shù)與粒子的均方根位移(MSD)成正比，表達(dá)式為[41]

其中，N是粒子數(shù)目，ri(t)-ri(0) 是第i個(gè)粒子在時(shí)間間隔t內(nèi)的運(yùn)動(dòng)位移.

2.2.4 氫鍵壽命

通常氫鍵采用純幾何準(zhǔn)則進(jìn)行定義:若兩個(gè)水分子的氧原子 (O···O) 間的距離小于0.35 nm，O···O—H 之間的夾角小于30°，則認(rèn)為兩水分子之間形成氫鍵.氫鍵的自相關(guān)函數(shù)為[42-44]

若一對(duì)水分子i和j在t時(shí)刻形成氫鍵，則h(t)的值為1，否則為0.氫鍵自相關(guān)函數(shù)C(t) 統(tǒng)計(jì)在0 時(shí)刻形成氫鍵的兩個(gè)水分子在t時(shí)刻仍然成氫鍵的概率.氫鍵壽命通過對(duì)自相關(guān)函數(shù)擬合得到，常用指數(shù)函數(shù)C(t)exp(-t/τ) 進(jìn)行擬合，式中τ為氫鍵的壽命[45].

3 結(jié)果和討論

理論計(jì)算結(jié)果表明，不同水模型計(jì)算得到的體相水的吸收譜與實(shí)驗(yàn)測(cè)量所得的吸收譜均存在一定程度的偏差，并且體相水的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收譜對(duì)溫度敏感 (如圖2 所示).首先，計(jì)算由不同水模型(SPC，SPC/E，TIP3P，TIP4P，TIP4P-Ew)構(gòu)成的體相水的太赫茲吸收譜，并將結(jié)果與實(shí)驗(yàn)測(cè)量結(jié)果[46]進(jìn)行對(duì)比.在定性分析上，理論計(jì)算結(jié)果與實(shí)驗(yàn)測(cè)量結(jié)果具有相同性，但是在定量分析上它們存在差異性.具體來(lái)說，如圖2(a)所示，由SPC，SPC/E，TIP3P，TIP4P，TIP4P-Ew 水模型構(gòu)成的體相水的太赫茲吸收譜都具有氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰、水分子內(nèi)部彎曲振動(dòng)的吸收峰和水分子內(nèi)部氫氧鍵拉伸振動(dòng)的吸收峰，這就意味著這些水模型能夠定性地描述體相水的各種振動(dòng)模式.定量上，雖然由SPC 和SPC/E 水模型構(gòu)成的體相水的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰的中心頻點(diǎn)在(20.11 ±0.06) THz 和(21.83 ± 0.04) THz，略大于實(shí)驗(yàn)測(cè)量值 (19.5 THz)，但是這兩種水模型的彎曲(45—50 THz)和拉伸振動(dòng)的吸收峰(90—105 THz)的位置確與實(shí)驗(yàn)測(cè)量值接近(圖2(a)上部).相反的是，由TIP3P，TIP4P，TIP4P-Ew 水模型構(gòu)成的體相水的彎曲(60—65 THz)和拉伸振動(dòng)的吸收峰(110—120 THz)的位置遠(yuǎn)偏離實(shí)驗(yàn)測(cè)量值，但是它們的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰的中心頻點(diǎn)((18.52 ±0.06) THz 對(duì)應(yīng)于TIP3P 水模型;(19.28 ± 0.05) THz對(duì)應(yīng)于TIP4P 水模型;(20.77 ± 0.05) THz 對(duì)應(yīng)于TIP4P-Ew 水模型)卻更接近實(shí)驗(yàn)測(cè)量值(圖2(a)下部).以上結(jié)果表明，SPC 和SPC/E 相比于TIP3P，TIP4P 和TIP4P-Ew 更能反映水的單分子振動(dòng)特性，然而TIP3P，TIP4P 和TIP4P-Ew 相比于SPC和SPC/E 更能反映水的氫鍵網(wǎng)絡(luò)的集體振動(dòng)特性.因?yàn)樗莫?dú)特性質(zhì)，例如較高的沸點(diǎn)和熱容量，來(lái)自于其復(fù)雜的氫鍵網(wǎng)絡(luò)，所以研究水的復(fù)雜氫鍵網(wǎng)絡(luò)結(jié)構(gòu)對(duì)于確定水的性質(zhì)至關(guān)重要.接下來(lái)，比較由不同水模型組成的體相水的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收譜(如圖2(b) 所示).理論模擬結(jié)果表明，相比于SPC/E，SPC 水模型的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰的中心頻率發(fā)生了紅移.這歸因于SPC 的氫氧原子電荷量的絕對(duì)值比SPC/E 小 (如表1 所列)，減小的電荷會(huì)影響偶極矩，從而降低介電常數(shù)，進(jìn)而影響分子間的相互作用[47]，最終導(dǎo)致SPC 的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰相比SPC/E發(fā)生了紅移.因?yàn)镾PC 和SPC/E 水模型具有相同的鍵長(zhǎng)rH—O、鍵角 ∠ HOH 參數(shù)，所以它們彎曲和拉伸振動(dòng)的吸收峰的頻率是一致的.同樣地，TIP3P，TIP4P 與TIP4P-Ew 的氫氧原子電荷量的絕對(duì)值的大小依次是TIP3P ＜ TIP4P ＜ TIP4P-Ew，隨著電荷的增大，氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰逐漸發(fā)生藍(lán)移，它們彎曲和拉伸振動(dòng)的吸收峰的頻率是一致的，歸因于它們也采用了相同的鍵長(zhǎng)rH—O、鍵角 ∠ HOH參數(shù).接下來(lái)，通過觀察不同水模型的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收譜之間的差異 (圖2(b))，發(fā)現(xiàn)它們的中心頻率依次是TIP3P ((18.52 ± 0.06) THz) ＜ TIP4P((19.28 ± 0.05) THz) ＜ SPC ((20.11 ± 0.06) THz)＜ TIP4P-Ew ((20.77 ± 0.05)THz) ＜ SPC/E ((21.83± 0.04) THz)，這預(yù)示著水的氫鍵網(wǎng)絡(luò)的振動(dòng)吸收峰的中心頻率與水的相互作用之間存在一定的關(guān)系.

圖2 不同水模型建模的體相水的太赫茲吸收譜和溫度對(duì)其太赫茲吸收譜的影響 (a) 不同水模型建模的體相水的太赫茲吸收光譜以及實(shí)驗(yàn)測(cè)量所得體相水的太赫茲吸收譜之間的對(duì)比;(b) 對(duì)于不同水模型構(gòu)成的體相水，其氫鍵網(wǎng)絡(luò)的太赫茲吸收譜;(c) 不同溫度下，SPC/E 水模型建模的體相水的太赫茲吸收譜;(d) 不同溫度下，對(duì)于SPC/E 水模型構(gòu)成的體相水，其氫鍵網(wǎng)絡(luò)的太赫茲吸收譜Fig.2.THz absorption spectra of bulk water modeled by different water models and effects of temperature on its spectra:(a) Comparison of THz absorption spectra of bulk water under different water models and its spectra from experimental measurement;(b) THz absorption spectra of hydrogen bond network of bulk water under different water models;(c) THz absorption spectra of bulk water modeled by SPC/E water model at different temperatures;(d) THz absorption spectra of hydrogen bond network of bulk water modeled by SPC/E water model at different temperatures.

因?yàn)闇囟葧?huì)改變氫鍵網(wǎng)絡(luò)的結(jié)構(gòu)，進(jìn)而影響氫鍵網(wǎng)絡(luò)的相互作用.較高的溫度更利于破壞氫鍵網(wǎng)絡(luò)，降低氫鍵的相互作用.為了研究氫鍵網(wǎng)絡(luò)的振動(dòng)吸收峰的中心頻率與水的相互作用之間的關(guān)系，選擇特定的水模型(SPC/E)，進(jìn)一步計(jì)算了體相水在不同溫度下的太赫茲吸收譜.從圖2(c)可知，溫度對(duì)水的氫鍵網(wǎng)絡(luò)的振動(dòng)吸收峰有顯著的影響，然而對(duì)水的彎曲和拉伸振動(dòng)的吸收峰的影響有限.特別地，隨著溫度從300 K 逐漸升高到370 K，水的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰的中心頻率逐步發(fā)生了紅移(如圖2(d) 所示)，也就是300 K 下為(21.83 ±0.04) THz，310 K 下為(21.73 ± 0.05) THz，330 K下為(21.39 ± 0.05) THz，350 K 下為(21.03 ±0.06) THz，370 K 下為(20.81 ± 0.07) THz.因此，可以初步得出結(jié)論:水具有越強(qiáng)的氫鍵，其氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰的中心頻率越大.

為了進(jìn)一步研究水的氫鍵網(wǎng)絡(luò)強(qiáng)弱與其振動(dòng)吸收峰的中心頻率之間的關(guān)系，分析了體相水的結(jié)構(gòu)和動(dòng)力學(xué)性質(zhì)(圖3).首先，計(jì)算了由不同水模型構(gòu)成的體相水的徑向分布函數(shù)g(r)，如圖3(a)所示.可以看出，隨著r的增加，g(r)函數(shù)值逐漸減小并收斂到1.0，這說明水分子表現(xiàn)出短程有序、長(zhǎng)程無(wú)序的微觀結(jié)構(gòu)特征，符合液體分子的排列規(guī)則.此外，由SPC/E，TIP4P 和TIP4P-Ew 構(gòu)成的體相水的g(r)具有兩個(gè)峰，跟實(shí)驗(yàn)結(jié)果一致[48]，但由SPC，TIP3P 構(gòu)成的體相水的g(r)只有1 個(gè)峰.特別地，由不同水模型構(gòu)成的體相水的第一水合層的密度是不一樣的，具體表現(xiàn)為g(r)的第1 個(gè)峰的峰值高度的差異，其峰值高度的大小關(guān)系依次為TIP3P (2.86) ＜ TIP4P (3.13) ＜ SPC (3.18) ＜TIP4P-Ew (3.43) ＜ SPC/E(3.51);然而，由不同水模型構(gòu)成的體相水的g(r)的第1 個(gè)峰值所處的位置(r值)具有相反的關(guān)系:TIP3P (0.274 nm) ＞TIP4P (0.273 nm)=SPC (0.273 nm) ＞ TIP4P-Ew(0.272 nm) ＞ SPC/E (0.271 nm).值得注意的是，g(r)的第1 個(gè)峰值高度值越大以及g(r)的第1 個(gè)峰值所處的r值越小，水的氫鍵網(wǎng)絡(luò)振動(dòng)的太赫茲吸收峰的中心頻率越大.這歸因于水分子與周圍的第一水合層內(nèi)的水分子形成氫鍵，g(r)的第1 個(gè)峰的峰值高度越大以及它的r越小，說明水分子與第一水合層內(nèi)的水分子結(jié)合得越緊密，從而具有越強(qiáng)的氫鍵相互作用，對(duì)周圍水分子的束縛也越強(qiáng).粒子的擴(kuò)散系數(shù)能夠說明周圍粒子對(duì)其的束縛情況.然后，計(jì)算了由不同水模型構(gòu)成的體相水的擴(kuò)散系數(shù)D，如圖3(b)所示，各D值的大小關(guān)系為TIP3P＞ TIP4P ≈ SPC ＞ TIP4P-Ew ＞ SPC/E，D的值依次是(3.97 ± 0.09) ×10—5，(2.36 ± 0.04) × 10—5，(2.34 ± 0.07) × 10—5，(1.25 ± 0.07) × 10—5和(1.12 ± 0.01) × 10—5cm2/s.將各D值的大小關(guān)系與其氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰的中心頻率的大小關(guān)系聯(lián)系起來(lái)，發(fā)現(xiàn)水的氫鍵對(duì)水的束縛越強(qiáng) (D值越小)，水的氫鍵網(wǎng)絡(luò)振動(dòng)的太赫茲吸收峰的中心頻率越大.同樣地，基于SPC/E 水模型，計(jì)算了溫度對(duì)體相水g(r)和D的影響.如圖3(c) 所示，隨著溫度的升高，g(r)的第1 個(gè)峰值高度不斷減小(300 K 下 為3.51，310 K 下 為3.39，330 K 下 為3.20，350 K 下為3.01，370 K 下為2.88)，它的第1個(gè)谷值卻不斷升高，以及第2 個(gè)峰的峰值也逐步降低甚至消失，這表明水的長(zhǎng)程有序性逐漸減小.不僅如此，隨著溫度升高，體相水的D逐漸增大，如圖3(d) 所示，表現(xiàn)為:D的值在300 K 時(shí)為(1.12 ±0.01) × 10—5cm2/s，310 K 時(shí)為 (1.54 ± 0.07) ×10—5cm2/s，330 K 時(shí)為 (2.47 ± 0.07) × 10—5cm2/s，350 K 時(shí)為 (3.54 ± 0.03) × 10—5cm2/s，370 K 時(shí)為 (4.80 ± 0.12) × 10—5cm2/s.這是因?yàn)闇囟壬呤沟盟肿拥倪\(yùn)動(dòng)加劇，削弱了水分子之間的氫鍵.以上結(jié)果表明:水分子具有越緊束縛的氫鍵網(wǎng)絡(luò)，其氫鍵網(wǎng)絡(luò)的振動(dòng)吸收峰的中心頻率越大.

圖3 不同水模型構(gòu)成的體相水的結(jié)構(gòu)和動(dòng)力學(xué)性質(zhì)的分析 (a) 不同水模型構(gòu)成的體相水的徑向分布函數(shù);(b) 不同水模型構(gòu)成的體相水的擴(kuò)散系數(shù);(c) 特定的水模型(SPC/E)下，溫度對(duì)體相水的徑向分布函數(shù)的影響;(d) 不同溫度下SPC/E 的擴(kuò)散系數(shù)Fig.3.Analyses of structure and dynamic behavior of bulk water modeled by different water models:(a) Radial distribution function of bulk water under different water models;(b) diffusion coefficient of bulk water under different water models;(c) under specific water model (SPC/E)，the effect of temperature on the radial distribution function of bulk water;(d) diffusion coefficient of bulk water at different temperatures.

最后，定量研究水的氫鍵強(qiáng)弱與其振動(dòng)吸收峰的中心頻率之間的關(guān)系(如圖4 所示).水分子之間的氫鍵處于動(dòng)態(tài)變化，它可能在某一時(shí)刻斷裂，隨后又與其他水分子形成氫鍵，氫鍵動(dòng)態(tài)變化導(dǎo)致體系總偶極距的變化，因此，根據(jù)漲落耗散定理[49]，氫鍵動(dòng)態(tài)變化的快慢將直接影響氫鍵網(wǎng)絡(luò)振動(dòng)的太赫茲吸收譜.氫鍵的平均壽命能夠衡量氫鍵動(dòng)態(tài)變化的快慢，它由氫鍵的自相關(guān)函數(shù)C(t)擬合得到.如圖4(a)所示，不同水模型的C(t)都近似于指數(shù)衰減，但是它們衰減的快慢不一樣，表明氫鍵壽命大小不一樣.從圖4(b)可知，SPC/E，TIP4P-Ew，SPC，TIP4P，TIP3P 的氫鍵壽命依次為(7.09 ±0.03)，(5.46 ± 0.02)，(3.23 ± 0.02)，(3.11 ± 0.02)和(1.66 ± 0.01) ps.將由不同水模型構(gòu)成的體相水的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收譜的中心頻率v與它們的氫鍵壽命τ聯(lián)系起來(lái) (如圖4(c)所示)，發(fā)現(xiàn)它們之間呈線性關(guān)系va1·τ+b1，其中a1=(0.37 ±0.05)，b1=(18.60 ± 0.22).不僅如此，氫鍵壽命隨著溫度升高逐漸減小，如圖4(d) 所示，體相水的τ在300 K 時(shí)為(7.09 ± 0.03) ps，310 K 時(shí)為(5.14 ±0.02) ps，330 K 時(shí)為(3.27 ± 0.01) ps，350 K 時(shí)為(2.32 ± 0.01) ps，以及370 K 時(shí)為(1.70 ± 0.01) ps.同樣地，將不同溫度下的體相水的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收譜的中心頻率v與它們的氫鍵壽命τ聯(lián)系起來(lái)(如圖4(e)所示)，發(fā)現(xiàn)τ越大，v也越大，它們之間也呈線性關(guān)系va2·τ+b2，其中a2=(0.33 ±0.04)，b2=(19.78 ± 0.12).值得注意的是，a1≈a2，這表明體相水的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收譜的中心頻率與它的氫鍵壽命之間的線性關(guān)系具有普適性.更進(jìn)一步剖析其中的機(jī)理，將水分子之間的氫鍵比作彈簧 (如圖4(f) 所示)，氫鍵網(wǎng)絡(luò)中的水分子之間的氫鍵相互作用越弱，說明彈簧的剛性越軟，即勁度系數(shù)k越小，具體表現(xiàn)為水分子對(duì)周圍水分子的束縛越弱，水分子的隨機(jī)擴(kuò)散越強(qiáng)，根據(jù)彈簧簡(jiǎn)諧振動(dòng)的頻率計(jì)算公式ω可知，ω隨k的減小而減小，表現(xiàn)為氫鍵網(wǎng)絡(luò)振動(dòng)的太赫茲吸收峰發(fā)生紅移.

圖4 氫鍵動(dòng)力學(xué)與氫鍵網(wǎng)絡(luò)振動(dòng)的太赫茲吸收峰的中心頻率之間的關(guān)系 (a) 由不同水模型構(gòu)成的體相水的氫鍵的自相關(guān)函數(shù);(b) 由不同水模型構(gòu)成的體相水的氫鍵壽命;(c) 不同水模型下，氫鍵網(wǎng)絡(luò)振動(dòng)的太赫茲吸收譜的中心頻率與氫鍵壽命的對(duì)應(yīng)關(guān)系;(d) 不同溫度下，體相水的氫鍵壽命;(e) 不同溫度下，水的氫鍵網(wǎng)絡(luò)振動(dòng)的太赫茲吸收峰的中心頻率與氫鍵壽命的關(guān)系;(f) 氫鍵網(wǎng)絡(luò)示意圖Fig.4.Relationship between hydrogen bond dynamics and the center frequency of THz absorption spectra for the vibration of hydrogen bond network:(a) Hydrogen bond autocorrelation functions of bulk water under different water models;(b) lifetime of hydrogen bond of bulk water under different water models;(c) for bulk water under different water models，the relationship between the center frequency of THz absorption spectra for the vibration of hydrogen bond network and the lifetime of hydrogen bond;(d) lifetime of hydrogen bond for bulk water at different temperatures;(e) at different temperatures，the relationship between the center frequency of THz absorption spectra for the vibration of the hydrogen bond network and the lifetime of the hydrogen bond;(f) schematic diagram of the hydrogen bond network.

4 結(jié)論

基于經(jīng)典分子動(dòng)力學(xué)模擬，首先研究了不同水模型的太赫茲吸收譜，以及溫度對(duì)太赫茲吸收光譜的影響.發(fā)現(xiàn)不同水模型均能夠定性地描述體相水的各種振動(dòng)模式，然而定量分析振動(dòng)吸收峰對(duì)應(yīng)的太赫茲頻率值卻有差異.進(jìn)一步研究了水的氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰的中心頻率與水的相互作用之間的關(guān)系，發(fā)現(xiàn)水的氫鍵相互作用越強(qiáng)，其氫鍵網(wǎng)絡(luò)的振動(dòng)吸收峰的中心頻率越大，即表現(xiàn)出藍(lán)移現(xiàn)象.為了更進(jìn)一步研究水的氫鍵網(wǎng)絡(luò)強(qiáng)弱與其振動(dòng)吸收峰的中心頻率的關(guān)系，計(jì)算了徑向分布函數(shù)和擴(kuò)散系數(shù)，研究結(jié)果表明水分子具有越緊束縛的氫鍵網(wǎng)絡(luò)，其氫鍵網(wǎng)絡(luò)的振動(dòng)吸收峰的中心頻率就越大.最后，基于以上分析，發(fā)現(xiàn)水的氫鍵壽命的大小與氫鍵網(wǎng)絡(luò)振動(dòng)的吸收峰的中心頻率之間存在線性關(guān)系.這些發(fā)現(xiàn)將有利于認(rèn)識(shí)水中復(fù)雜的氫鍵網(wǎng)絡(luò)動(dòng)力學(xué)性質(zhì)和加深在分子尺度對(duì)水的振動(dòng)運(yùn)動(dòng)的理解，同時(shí)促進(jìn)太赫茲技術(shù)在生物醫(yī)學(xué)上的應(yīng)用.