2010年12月6日 星期一

Molecular Dynamic Simulation - 威佐、雙魚

--------- 2010 . 12. 06 ----------
威佐 :

目前已經將幾個常用的 Operations 取出,
需要計算的力包括:
Bonded(Bond length、Bond angle、Dihedral 三種力)和 Non-Bonded,

我相信這是簡化過的形式,但對於我們做加速運算已經足夠,

我現在使用的 Code 是學長給的那本書中的 Code (老師我沒有自己寫一個...)。
但是因為要先嘗試做 translator,所以先暫用書中的 Code,等到評估期過後,
我們可以自己改一份 MD 的 Code,

另外,現在成熟的 MD 套件大概是 CHAMBER、AMBER、NAMD、GROMACS,
可以下載使用,看一下別人如何做的,當然我們不需要做到像這些軟體這麼複雜。


雙魚:

著手將 MD 的 Code 改成 Multi-Thread,
目前先規劃將 Non-Bonded Force 的部份給平行,
如果修改完成,可以再加入 Bonded 的部分。

下一步,就是將它寫成 SPCI 的 Component。

2010年12月5日 星期日

ROC - 雙魚

=== 2010.11.23 ===

更改成matlab可連結的函式後
把此份.c檔放到要執行的matlab資料夾內

編譯.c檔的流程如下:
(1)如果是第一次用mex,先設定編譯器
執行 mex , 畫面會列出你電腦上現有的編譯器
選擇好後電腦會記下,以後使用就不用執行這步驟

(2)編譯.c檔 : 執行 mex xx.c (xx為檔名)
如果順利編譯通過,同一個資料夾下會出現xx.dll檔案
(3)接著就可以執行matlab code


比較原本matlab程式和呼叫 C library的matlab程式時間
(1)原本matlab程式:
for j = 1:1000 呼叫ROC計算函式 end
時間為: 8分鐘

(2)呼叫 C library的matlab程式:
for j = 1:1000 呼叫 ROC_scream() end
時間為: 0.679 s

呼叫 C library的matlab程式時間節省不少
11/23把此程式交給龔老師,也跟他說明了使用方法
ROC改寫到底為止

=== 2010.11.19 ===

把改寫後的C code包成matlab可用的library型式
即可在matlab程式中直接使用改寫後的code
參考連結: http://twntwn.info/blog/ajer001/archives/790
按照連結中的方法來更動原本的C code

matlab迴圈呼叫的寫法:
for j = 1:1000
[auc_resp_sf fa_resp_sf hr_resp_sf] = ROC_scream(all4type, all4respyesno_sf);
end

說明: ROC_scream()為改寫後的C code
all4type,all4respyesno_sf為輸入的數列
[auc_resp_sf fa_resp_sf hr_resp_sf]為output
fa_resp_sf和 hr_resp_sf是數列,auc_resp_sf為1*1的double

更動步驟如下:
(1)在code中開頭加入 #include "mex.h"
(2)將程式主要部分放在 mexFunction()中
mexFunction(int nlhs, mxArray *plhs[],
int nrhs, const mxArray *prhs[])
mexFunction()為程式的進入點,prototype不可隨意更動


(3)要把matlab輸入數列換成 C陣列
value = mxGetPr(prhs[0]);
sign = mxGetPr(prhs[1]);
bufsize = mxGetM(prhs[0]);
n = mxGetN(prhs[0]);

說明:
plhs: ROC_scream函式等號左邊參數本體
prhs: ROC_scream函式等號右邊參數本體
mxGetPr:取得矩陣之實數部分
mxGetM:取得矩陣之第一維大小

value會取得ROC_scream()中 all4type數列
sign會取得ROC_scream()中all4respyesno_sf數列
bufsize是all4type數列的row數量
n是all4type數列的column數量


(4)給予輸出陣列的記憶體空間
plhs[0]= mxCreateDoubleMatrix(mA,nA,mxREAL);
area = mxGetPr(plhs[0]);

說明:
利用mxCreateDoubleMatrix()建立空間
mA是row數量 nA是column數
plhs[0]代表為回傳的第一個參數,並取出取得矩陣之實數部分


=== 2010.11.08 ===

跟龔老師拿了他們使用輸入輸出的數據
以便於我們做更多數據的比對
每次輸入的a和b數列各160筆,要跑1000次

(1)和matlab code輸入相同的a,b數列
比對C code算出的面積數字是否和原本一樣
比對的結果皆正確

(2)比較兩份code時間 (迴圈1000次,每次160筆data)
matlab迴圈寫法為:
for j = 1:1000 呼叫ROC計算函式 end

時間結果如下:
Matlab時間為= 483.58 s ≈ 8分鐘
C code時間為= 1.984 s ≈ 約兩秒


=== 2010.10.30 ===

完成 ROC C code改寫

ROC程式內容如下:

(1)從輸入的數列中找出沒有重複的uv數列
input: 數列a和b (random ,各200個數值)
找出uv數列:為a和b中沒有重複的數值

舉例說明:
a = {-2,-1,0,3,8 }
b = {-3,1,3,0,8 }
uv = {-3,-2,-1,0,1,3,8}

(2)找 hit rate (hr) and false alarms (fa)
fa: a中有幾個大於等於uv數值 / a數列總合
hr: b中有幾個大於等於uv數值 / b數列總合
uv中每個數值都要檢查一次

舉例說明: uv ={ -3,-2,-1,0,1,3,8 }
a ={ -2,-1,0,3,8 } , b ={ -3,1,3,0,8 }
當uv[0]= -3 , fa[0]=1 and hr[0]=1
當uv[1]= -2 , fa[1]=1 and hr[1]=4/5


(3)計算ROC曲線下的面積 (參考圖)


matlab code : auc = trapz(fa,hr);
trapezoid function使用梯形公式算曲線下的面積
梯形的高為 fa[i] – fa[i+1]
上底為hr[i],下底為hr[i+1]
Auc = 梯形公式 ((上底+下底)*高)/2

舉例說明:
trapz([123],[124])
梯形第一塊: ((1+2)*(2-1))/2
梯形第二塊: ((2+4)*(3-2))/2


=== 2010.10.19 ===

在分析fMRI數據時,需要花費大量的計算和時間
心理系龔老師希望我們來幫忙減少code的計算時間
龔老師提供了他們主要計算部分的matlab程式碼
其中以"ROC"此部分所花費的時間最多

因此我們從此部分著手
老師的想法是把ROC此部分改成C code
改完後測時間看是否還要再作平行的方法