目的：在NVE集合中模擬泛素，並分析溫度分佈。
1、在終端窗口中，轉到目錄2-4-temp /。
對於本練習，您需要一個足夠長的模擬以很好地採樣足夠數量感興趣的樣本。在這裡，我們需要採樣動能以及溫度的波動。建議您不要自己執行這個模擬，你可以使用已提供的文件example-output/ubq-nve.log。但是，如果你執意選擇運行自己的模擬，則可以在下面找到所有必要的文件和說明。即使不執行模擬，我們也建議您閱讀整個部分。
BOND：232.561128
ANGLE：690.805988
DIHED：397.649676
IMPRP：41.250136
ELECT：23674.863868
VDW：1541.108344
BOUNDARY：0.0
MISC：0.0
KINETIC：4579.147764
TOTAL：16192.340832
TEMP：313.429264
TOTAL2：16173.821268
TOTAL3：16174.663916
TEMPAVG：313.140436
PRESSURE：-47.513352
GPRESSURE：-19.105588
VOLUME：71087.006196
PRESSAVG：3.340104
GPRESSAVG：3.4207
溫度：一個係綜的溫度由一下確定：
這個表達式來自於應用於動能的統計力學的均分定理
N是系統的原子數量
該模擬在微正則係綜（NVE，即常數，體積和能量）下執行。為啟動NAMD運行而提供的配置文件稱為ubq-nve.conf。該配置文件的主要功能是：
①該模擬從單元1中執行的水分子球平衡模擬中獲取重新restart文件。這些文件位於目錄1-2-sphere /中，並命名為ubq_ws，即ubq_ws.restart.vel等。
②初始溫度由速度restart文件確定，如上式所示。此初始溫度對應於Ť〜313 K，是水球體模擬的溫度。
③時間步長為2 fs；如果打開了rigidBonds，則此時間步是可以接受的。
④仿真將運行1 ns。
2運行分子動力學模擬
namd2 ubq-nve.conf > ubq-nve.log &
3、如果您自己沒有運行模擬，則需要從提供的log日誌文件中獲取溫度數據。從目錄example-output /複製此日誌文件。
copy example-output\ubq-nve.log。
如果您自己運行模擬，則日誌文件將由NAMD生成，並且您無需複制它。
4、在VMD TkCon窗口中，通過鍵入cd ../2-4-temp/導航至2-4-temp /目錄。
5、為了從日誌文件中獲取溫度數據，我們將再次使用腳本namdstats.tcl，該腳本已由我們提供。通過鍵入以下內容：
data_time TEMP ubq-nve.log
現在，每個時間步長和系統的相應溫度都在文件TEMP.dat中。該腳本將包含所有時間步，可能需要一些時間才能完成。請注意，系統沒有最小化，因為這已在初始模擬中執行。
為了繪製數據圖，可以使用Excel。
6、打開Excel。
7、單擊文件→ 打開...，導航到目錄2-4-temp /，然後選擇文件TEMP.dat。確保下拉菜單文件類型設置為“所有文件”。
8、出現“文本導入嚮導”時，為文件類型選擇“分隔符”，然後單擊“下一步”。在“分隔符”類別下，單擊“Tab鍵”。 “數據預覽”窗口應將您的時間步長和溫度數據分為兩列。單擊完成。時間步長和溫度值將分別出現在列A和B中。
我們將構建溫度值的直方圖，並對其擬合高斯分佈。
9、構造直方圖的一列柱標籤。在單元格C1中，輸入數字300，在單元格C2中，輸入公式“= C1 + 0.25”，然後按Enter鍵。單擊單元格C2，然後單擊並按住其右下角，然後向下拖動到第121行並釋放。您應該在C列中有300至330的溫度值列表。這是我們的直方圖分區。
10、現在，我們將使用“數據分析”工具來創建直方圖數據。點擊數據→ 數據分析...。如果該選項不存在，則跟我之前說的方法一樣。並從列表中選擇分析工具庫來添加數據分析工具。
11、從數據分析窗口中選擇直方圖，然後單擊確定。輸入以下單元格範圍，然後單擊“確定”：
輸入範圍：$B$1：$B$10000
接收範圍：$C$1：$C$121
將使用您的溫度柱和給定溫度柱中發生的溫度測量次數創建一個新表。我們將繪製這些值。
12、通過單擊插入→ 圖表...，選擇柱狀圖作為圖表類型，選擇左上角的第一個圖表子類型，然後單擊“下一步”，來繪製數據的直方圖。
12、通過單擊插入→ 圖表...，選擇柱狀圖作為圖表類型，選擇左上角的第一個圖表子類型，然後單擊“下一步”，來繪製數據的直方圖。
13、在步驟2中，單擊“選擇數據”選項卡。刪除第一個系列，然後單擊“添加(A)”按鈕添加一個新系列。單擊“系列值：”字段中的箭頭按鈕，選擇列B中的每個數據值（行2-122），然後單擊確定完成操作。對“水平（分類）（X）軸標籤執”行相同的操作：使用A列的2-122行。單擊“完成”。
14、隨意修改圖表。我們建議右鍵單擊直方圖柱，選擇“數據系列格式”，單擊“選項”選項卡，並將“間隙寬度：”設置為0。
這種分佈看起來是否熟悉？您的分佈看起來應該像高斯分佈。
溫度波動：動能En以不同的形式進行麥克斯分佈
由以下推導可得：
總動能Ek = Pj 1 2 mjvj2的分佈，根據中心極限定理，近似高斯分佈。
溫度 T = 2Ek/3kB的分佈函數，波動是∆T = T - T0
注意，對於N趨近於無窮大時，分佈是很陡峭的，但是對於一個有限的系統來說是很平坦的。對於T=300K和N=10001有一個σ-8K。
15、現在，您將用正態分佈擬合該分佈。我們將以類似於Maxwell-Boltzmann能量分佈的方式進行此操作。
16、在單元格C2中，輸入公式：
= D$2 * EXP（-（（A2-D$3）^2 ）/ D$4）
D列中的值是高斯分佈的參數。在單元格D2，D3和D4中分別輸入初始猜測400、315和100。然後將C2中的公式向下擴展到單元格C122。這將創建高斯分佈來擬合原始數據。該參數對應於歸一化常數，平均溫度，和σ 2。
17、在柱狀圖，在柱狀圖單擊鼠標右鍵，並點擊“選擇數據”。單擊“添加（A）”選項卡，然後在列C的第2-122行中添加帶有“值”的新系列。單擊確定。
18、右鍵單擊剛創建的紅色直方圖條，選擇更改圖標類型，選擇折線圖作為圖標類型，然後選擇左上角的第一個Chart子類型。單擊確定。
您可能想要更改高斯擬合的線寬，可以通過右鍵單擊直方圖中的線並選擇“格式化數據序列...”來進行。
我們繪製了高斯曲線，但擬合度不是最佳的。如果更改參數列D的任何值，我們將看到擬合改變。現在，我們將確定最佳參數。
19、在單元格E2中，鍵入公式=（C2-B2）^ 2.然後單擊並按住單元格E2的右下角和向下拖動到行122填補對每個段的誤差值。
20、在單元格F2中鍵入= sum（E2：E122）。這是我們擬合的總平方誤差，也是我們尋求使用規劃求解工具使之最小化的原因。
21、單擊數據→ 規劃求解。如果該選項不存在，則需要通過單擊“工具” →“ 添加項”並從列表中選擇“求解器加載項” 來添加“求解器”工具。
22、在規劃求解參數窗口中，設置以下單元格範圍，然後單擊
解決：
設置目標單元格：$F$2
等於：最小值
通過更改單元格：$D$2：$D$4
Solver將找到最適合我們原始溫度數據的高斯分佈參數值。
23、將通過此方法獲得的溫度平均值與通過在VMD TkCon窗口中鍵入namdstats獲得的平均值進行比較：data_avg ubq-nve.log
24、現在，看看偏差σ 2 = 2 Ť 2 / 3 Ñ，可以從你的數據來計算。注意它如何取決於系統的大小！
提示：您可能會從該示例中註意到，一方面，單個蛋白質的行為類似於無限的熱力學係綜，從而再現各自的平均溫度。另一方面，它們以溫度波動的形式顯示出其有限性的跡象。通常，有趣的是，原子的速度以及蛋白質的動能充當其溫度計！
25、如果需要，保存在Excel中創建的文件，然後關閉Excel。
沒有留言:
張貼留言