以多圖形處理器加速近似字串比對

全文

(1)國立臺灣師範大學科技應用與人力資源發展學系碩士論文. 以多圖形處理器加速近似字串比對. 研究生：張浩平指導教授：林政宏. 中華民國一 ○ 一年七月.

(2)

(3) 謝誌研究生活終於要告一個段落了，回想起這三年的點點滴滴，有回憶不完的酸甜苦辣。不過，完成學業的喜悅也正式邁向另一階段的動力，接下來路途無論如何，只要堅持下去，一定可以得到成果。一路走來，有許許多多人的支持與陪伴，點滴在心頭。感謝家人的支持，讓我可以專心的完成學業；感謝朋友與同學的陪伴，在遇到困境的時候，可以幫助我紓解情緒；感謝系上老師的教導，讓我得到更多方面的知識；感謝系辦的倩綾、曉鈴、秋逸與潔懋四位助教，妳們幫了我很多忙，真的很感謝妳們的容忍與體諒。還有HCLab的學弟們，不管是程式還是伺服器系統上的問題，有你們可以互相討論指教，我才可以解決問題，感謝你們，也祝你們順利。感謝蕭顯勝老師以及陳勇志老師，在我緊張到不行的口試中，用幽默的對話方式化解我緊張的情緒；最後，感謝林政宏老師，提供我非常多的資源與協助，讓我可以找出研究的方向，順順利利的畢業。也因為老師熱心的指導，讓學生有更多的機會去學習不同領域的知識。希望大家都可以事事順利，天天開心。張浩平於屏東中華民國一○一年七月. i.

(4) ii.

(5) 以多圖形處理器加速近似字串比對研究生：張浩平指導教授：林政宏. 中文摘要字串比對在許多領域中被廣泛運用，例如網路入侵偵測系統、網際網路搜尋、去氧核糖核酸序列比對等等；其中，字串比對可區分為固定字串比對與近似字串比對兩類。固定字串比對是指找出所有字串樣式於輸入文字中出現的位置，所有字串樣式必須精確的比對，不容許任何錯誤；而所謂近似字串比對是指所搜尋的字串樣式則可經由插入、刪除及替代等有限次數的動作，轉換成輸入文字中的某部分。近似字串比對的演算法可區分為 dynamic programming 與 bit-parallelism ； Dynamic programming需經龐大運算及記憶體空間記錄誤差值，故在處理大量資料時，將為此演算法之瓶頸；反之，bit-parallelism運用邏輯運算子模擬非確定有限狀態自動機進行比對，速度快且節省記憶體。近似字串比對的運算量大且非常耗費時間，尤其在針對大量的輸入文字比對大量的字串樣式時，其耗費的時間更為明顯。本研究將分析並實現bit-parallelism與dynamic programming於NVIDIA GPU上，實驗結果顯示在處理2Gbytes的輸入文字時，執行於單一個GPU的bit-parallelism較執行於單一執行緒CPU版本的bit-parallelism快上7倍的加速。本研究並進一步透過openMP連結多個GPU的加速，其結果顯示在處理2Gbytes的輸入文字時，以2個GPU加速之bit-parallelism較執行於單一執行緒CPU版本的bit-parallelism快上10倍的速度。關鍵詞：近似字串比對、GPU、dynamic programming、bit-parallelism iii.

(6) iv.

(7) Accelerating Approximate String Matching Using GPUs Author: Chang, Hao-Ping Adviser: Lin, Cheng-Hung. ABSTRACT Sting matching has been widely used in many areas, such as NIDS, web searching, and DNA sequence matching, etc. Sting matching can be classified into two main categories, exact string matching and approximate string matching. Exact string matching finds all occurrences of given patterns in an input string while approximate string matching allows difference between a string and a pattern caused by insertion, deletion and substitution. Approximate string matching has two main approaches, dynamic programming and bit-parallelism. Because dynamic programming needs very large memory space for store results, it is difficult to deal with huge inputs and not proper to be implemented on GPUs. On the other hand, as bit-parallelism uses bitwise operators to simulate NFA, it is fast and memory-efficient. Because approximate string matching is a very compute-intensive task, accelerating approximate string matching has become critical for processing big data. In this thesis, dynamic programming and bit-parallelism are implemented and evaluated on NVIDIA GPUs. The experimental results show that the bit-parallelism performed on GPUs achieves 7 times faster than single threaded CPU version for processing 2Gbytes data. Furthermore, multiple GPUs are adopted to increase throughput of processing big data. For processing 2Gbytes data, adopting two GPUs can achieve 10 times faster than single threaded CPU version. Keywords: approximate string matching, GPU, dynamic programming, bit-parallelism v.

(8) vi.

(9) 目錄謝誌 .................................................................................................................. i 中文摘要 ........................................................................................................ iii ABSTRACT .................................................................................................... v 目錄 ............................................................................................................... vii 表次 ................................................................................................................ xi 圖次 .............................................................................................................. xiii 第一章. 緒論 .............................................................................................. 1. 第一節研究背景 .................................................................................. 1 第二節研究動機與目的 ...................................................................... 3 第三節文獻探討 .................................................................................. 3 第四節研究流程 .................................................................................. 5 第五節研究貢獻 .................................................................................. 5 第二章. 近似字串比對演算法 .................................................................. 7. 第一節 Dynamic programming ............................................................. 7 第二節 Bit-parallelism ........................................................................ 10 第三節綜合分析 ................................................................................ 15 vii.

(10) 第三章. 平行處理架構 ............................................................................ 17. 第一節任務並行 ................................................................................ 17 第二節邊界偵測的問題 .................................................................... 18 第三節數據並行 I .............................................................................. 19 第四節數據並行 II............................................................................. 20 第五節本章小結 ................................................................................ 21 第四章. 實驗結果 .................................................................................... 23. 第一節實驗環境 ................................................................................ 23 第二節以 CPU 加速........................................................................... 24 第三節以 GPU 加速 .......................................................................... 26 第四節以多個 GPU 加速 .................................................................. 27 第五章. 結論 ............................................................................................ 31. 參考文獻 ....................................................................................................... 33 附錄 ........................................................................................................... 37 附錄一、Dynamic programming（CPU 版本） .................................. 39 附錄二、Bit-parallelism（CPU 版本） ............................................... 44 附錄三、Bit-parallelism（GPU 版本） ............................................... 50 viii.

(11) 附錄四、Bit-parallelism（Multi-GPU 版本） ..................................... 58. ix.

(12) x.

(13) 表次表 2-1 Dynamic programming 與 bit-parallelism 比較表 ........................... 16 表 4-1 Host 與 device 主機設備規格表 ..................................................... 23 表 4-2 CPU 加速 dynamic programming 之數據 ....................................... 24 表 4-3 CPU 加速 bit-parallelism 之數據 .................................................... 24 表 4-4 GPU 加速 bit-parallelism 之數據 ................................................... 26 表 4-5 Multi-GPU 加速 bit-parallelism 之數據.......................................... 28. xi.

(14) xii.

(15) 圖次圖 2-1 初始化 dynamic programming 以二維表格表示 ............................. 8 圖 2-2 計算 dynamic programming ............................................................. 9 圖 2-3 透過 back-traverse 追蹤過程 ......................................................... 10 圖 2-4 非確定有限狀態自動機（NFA） ................................................. 11 圖 2-5 說明如何儲存𝐷0 與𝐷1 ................................................................. 12 圖 2-6 設定𝑀𝑇𝑗遮罩 ................................................................................. 14 圖 2-7 模擬 bit-parallelism 內部 NFA 運作流程 ...................................... 15 圖 3-1 任務並行架構 ................................................................................... 18 圖 3-2 邊界偵測問題 ................................................................................ 18 圖 3-3 防止邊界偵測問題 ........................................................................ 19 圖 3-4 數據並行 I...................................................................................... 20 圖 3-5 數據並行 II .................................................................................... 21 圖 4-1 以 CPU 加速 dynamic programming 效能直條圖 .......................... 24 圖 4-2 以 CPU 加速 bit-parallelism 效能直條圖........................................ 25 圖 4-3 以 GPU 加速 bit-parallelism 效能直條圖 ....................................... 27 圖 4-4 分段讀取檔案過大的輸入文字 ....................................................... 28 圖 4-5 實驗結果折線圖 ............................................................................... 29 xiii.

(16) xiv.

(17) 第一章. 緒論. 本章首先介紹固定字串比對與近似字串比對的差異，爾後收集分析相關論文，提出研究動機及研究流程。本章主要分為研究背景、研究動機、文獻探討、研究流程及研究結果共五小節。說明如下：. 第一節. 研究背景. 字串比對（string matching）在許多領域中被廣泛運用，例如網路入侵偵測系統（network intrusion detection system）、網際網路搜尋（web searching）、去氧核糖核酸（DNA）序列比對等等。而字串比對演算法可區分為固定字串比對（ exact string matching ）與近似字串比對（approximate string matching）兩類。其中，固定字串比對必須找出所有字串樣式（string pattern）於輸入文字（input stream）中出現的位置，所有字串樣式必須精確的比對，不容許任何錯誤。另一方面，近似字串比對所搜尋的字串樣式則可容許特定數量的誤差發生，其容許誤差就是允許比對的字串樣式可以經由插入（insertion）、刪除（deletion）及替代（substitution）等動作，轉換成輸入文字中的某一部分（Navarro, 2001）。這些動作有一定的容許誤差，一般稱為 edit distance 或 Levenshtein distance，此容許誤差由使用者決定。舉例來說，給定一個字串樣式為"survey"與一個輸入文字為"surgery"，以固定字串比對而言，"survey"沒有 match "surgery"，但以近似字串比對而言，"survey"可以經由替代 v 為 g，與插入一個 r 於 ey 之間，而 match "surgery"，一般稱此兩字串之間的 edit distance 為 2；換言之，在 edit distance 為 2 的條件下，"survey" matches "surgery"。由於容許特定數量的誤差，近似字串比對更廣泛被運用於網際網路 1.

(18) 搜尋（Qin, Li, Zhu, & Tian, 2010; Wu, Zuo, Hu, Zhu, & Li, 2008）、去氧核糖核酸序列比對（Ivanko, 2006; Li, Ni, Wong, & Leung, 2011）、拼字檢查（Lyras, Sgarbas, & Fakotakis, 2007; Soleh & Purwarianti, 2011）、漢明碼搜尋音樂（Dong & Qi, 2010; Liu, Huang, Yang, & Zhang, 2009）、封包過濾（Liu, Lin, Li, & Lee, 2005; Shi & Xie, 2011）等應用。而近似字串比對演算法可以大致區分為 dynamic programming（Kaplan & Kaplan, 2004; Smith & Waterman, 1981; Si, Yang, Lu, Sun, & Mei, 2009）與 bit-parallelism （Baeza-Yates & Navarro, 1999; Hyyrö, 2008; Myers, 1999; Wu & Manber, 1992）兩大類；其中 dynamic programming 是採用建表的方式來計算 edit distance，並以 back-traverse table 的方式來求得比對的位置，傳統的 dynamic programming 演算法時間複雜度與空間複雜度分別為 O（𝑚𝑛） + O（𝑚 + 𝑛）與 O（𝑚𝑛），其中 m 與 n 分別代表字串樣式的長度與輸入文字的長度；故記憶體需求便成為傳統 dynamic programming 在處理大量資料時的瓶頸。另一方面，bit-parallelism 利用 AND、OR、XOR 等邏輯運算模擬非確定有限狀態自動機（Nondeterministic Finite Automaton, NFA）的運作，比對是否有相似的字串樣式存在於輸入文字中，其時間複雜度與空間複雜度分別為 O（n）與 O（（m-k）（k+2）），k 為 edit distance（Baeza-Yates & Navarro, 1999）。Bit-parallelism 因為採用大量的位元運算子（bitwise operation），因此效率比 dynamic programming 要高，記憶體需求亦遠小於 dynamic programming 演算法，唯受限於長度較短的字串樣式。上述兩種演算法，有其各別適用的時機，例如 dynamic programming 可提供在不同比重下的近似字串比對，像是 smith-waterman 演算法除了利用 dynamic programming 的核心運算外，還須搭配核酸與蛋白質序列比對工具（Basic Local Alignment Search Tool, BLAST）加重不同元素的 2.

(19) 比重，而 bit-parallelism 便不適用於此項運算。Bit-parallelism 的優勢在於處理速度快與節省記憶體，因此在需要快速搜尋字串樣式時， bit-parallelism 便可以提供較 dynamic programming 省時的方法。. 第二節. 研究動機與目的. 近似字串比對的運算量大且非常耗費時間，尤其在針對大量的輸入文字比對大量的字串樣式時，其耗費的時間更為明顯。在近似字串比對演算法中，由於 dynamic programming 必須透過建表的方式來進行運算，其表格大小為整個輸入文字與字串樣式相乘的二維陣列；因此在處理一個 1Gbytes 的輸入文字與長度為 10 的字串樣式時，所佔據的記憶體將近 10Gbytes，一般個人電腦是無法處理的。反之，bit-parallelism 利用邏輯運算來模擬非確定有限狀態自動機，在實現此演算法時，除了將輸入文字搬入硬碟之外，只需儲存非確定有限狀態自動機，此狀態自動機利用一個無號長整數型態來表示，大小為 4 個位元組（Bytes）。而非確定有限狀態自動機使用大量的位元向量表示與運算，也較貼近 CPU 與 GPU 的核心運作模式，除了所需的記憶體也不大之外，其處理速度也較快。綜合以上因素，為了降低時間與空間成本，本研究的目的為探討如何運用多圖形處理器（Graphical Processing Units, GPU）加速近似字串比對之演算法。. 第三節. 文獻探討. 近似字串比對應用範疇種類繁多，從比對字母到比對圖形、影像皆有其應用。由於容許誤差值的特性，可彌補資訊不充足的情況下，搜尋到期待的資料。. 3.

(20) 生物資訊方面，去氧核醣核酸序列比對可以得知兩基因序列的相似度，可幫助建置生物序列資料庫；在取得一個未知物種之基因序列片段時，可透過去氧核醣核酸序列比對在生物序列資料庫中，搜尋出在目前所知的物種中，相似度最高的生物。基因檢測為另一近似字串比對之應用，可提供一個檢測先天遺傳疾病或癌症的檢測（Dinu & Ionescu, 2011）。對於網際網路的應用，近似字串比對提供搜尋引擎更強大的搜尋能力，方面使用者得到多種不同組合的資訊。資訊安全方面，傳統的網路入侵偵測是採用固定字串比對，但攻擊行為往往夾雜許多雜訊，所以近似字串比對可以當作過濾器，來偵測隱藏的攻擊行為（Song, Xue, & Wang, 2006）。而 GPU 多執行緒的技術提供強大的運算能力，可以更快速的獲取想要的結果與資訊。許多研究便透過 GPU 加速 smith-waterman 演算法， smith-waterman 是由 dynamic programming 衍伸而來的演算法（Liu, Schmidt, Voss, Schroder, & Muller-Wittig, 2006; Liu, Huang, Johnson, & Vaidya, 2006; Munekawa, Ino, & Hagihara, 2008）。而在立體影像的領域，也有許多利用 GPU 加速基於 dynamic programming 來實現影像的演算法（Moslah, Valles-Such, Guitteny, Couvet, & Philipp-Foliguet, 2009; Kalarot, Morris, & Gimel’farb, 2010; Liao, Lin, & Medioni, 2011）。然而本研究發現以 bit-parallelism 為加速對象之研究相對較少。另外 bit-parallelism 透過 AND、OR、XOR 等邏輯運算，比起 dynamic programming 需要作大量的加法與條件判斷式來說，較符合 GPU 等硬體機器的運作。結合以上種種因素，本研究欲以 bit-parallelism 為研究對象，透過多個 GPU 來加速效能。. 4.

(21) 第四節. 研究流程. 本研究流程主要分為三個階段，首先本研究分析 dynamic programming 與 bit-parallelism 兩種演算法，瞭解其差異與適用時機，並分析此兩種演算法平行加速的可行性。接下來，本研究將各別實現此兩種演算法於 CPU 與 GPU 上，並比較 CPU 與 GPU 的成效。最後本研究將給予一個大量的輸入文字，透過多個 GPU 加速 bit-parallelism。. 第五節. 研究貢獻. 本研究初步實驗於四核心 CPU 主機並結合 pthread library 實現了 dynamic programming 與 bit-parallelism 兩種近似字串比對主要的演算法，以輸入文字 256Mbytes 為例，實驗結果分別都比傳統單一執行緒 CPU 版本的 dynamic programming 及 bit-parallelism 快上 4 倍與 5 倍左右的速度；另外在處理 2Gbytes 的輸入文字時，以單一個 GPU 加速的 bit-parallelism 可得到比傳統使用單一執行緒 CPU 版本的 bit-parallelism 快上 7 倍左右的加速，而比以 pthread library 加速的 bit-parallelism 快上 2 倍的速度。而本研究最後將進一步利用多個 GPU，結合 openMP 技術提升近似字串比對演算法的速度。. 5.

(22) 6.

(23) 第二章. 近似字串比對演算法. 前一章簡述了字串比對的種類，近似字串比對與固定字串比對最大的差異就是容許誤差。而實務應用中，使用者並非擁有詳盡的資訊，以網際網路搜尋引擎為例，輸入關鍵字師範大學的結果，可能有台灣師範大學、屏東師範大學或彰化師範大學等等，其應用較為廣泛。因此，本章將著重於近似字串比對相關的演算法簡介，並配合範例來解說其運作模式。本章分別介紹 dynamic programming、bit-parallelism 兩種演算法與綜合分析，茲說明如下：. 第一節. Dynamic programming. Dynamic programming是非常著名的近似字串比對演算法，其運用範圍非常廣泛，尤其dynamic programming可以給予使用者非常詳細的資訊，往後也有許多演算法是基於dynamic programming發展而來的。傳統dynamic programming主要分為兩個步驟，第一個步驟為找出 edit distance，另一步驟為back-traverse table。一. Edit distance i. ii.. 𝑗. iii.. { 1 + 𝑚 𝑛(. (. ). ). 利用上述演算式得以計算edit distance，首先要求得字串樣式與輸入文字兩個字串的edit distance，透過上述的方程式，將字串樣式以表示，輸入文字以表示，以建置一個二維向量的表格。舉例來說，假設有一個字串樣式為kittern，另有一個輸入文字為sitting，則利用運算式i與ii做 7.

(24) 初始化的動作。如圖2-1所示：. 圖 2-1 初始化 dynamic programming 以二維表格表示. 初始化表格後，接著要循序漸進地透過運算式iii之結果來填滿表格。在. 的位置時，首先判別字串樣式的第i個字元是否與輸入字串的第j個. 字元相同，若相同則將左上方的值填入此位置；若不相同則將正上方、左上方位置。以. 、左前方. 三處的值，取其最小值加1後填入此. 為例，因為為與為不相同，且左上方的值為最小值，. 因此將其加1後填入此位置，最後得. 1。其他位置則以此類推，但. 需依序由左至右，由上往下填表，索引值參照C程式語言從0開始。. 8.

(25) 圖 2-2 計算 dynamic programming. 二. Back-traverse table 所謂back-traverse table即是利用圖2-2之結果，往回追蹤並記錄字串樣式轉換成輸入文字之過程。經由前一步驟填滿表格後可以得到. ，. 換句話說，字串樣式kitten可透過三次前述之插入、刪除、替代等動作轉換成輸入文字sitting，其追蹤準則如iv所示。 iv. 在運算式iv中，. 𝑚 𝑛(. ). 表示最右下角之位置，即是指圖2-3中. 之位置，. 而左前方位置的值為最小值，故儲存其走向並繼續計算走向直到止。. 9. 為.

(26) 圖 2-3 透過 back-traverse 追蹤過程. 圖2-3箭頭所指方向即是所儲存的走向，倘若往左前方走，表示插入；往正上方走則表示刪除；若往左上方走，則須判別其值是否有減1，有則表示替代，沒有則表示match。因此圖2-3所呈現的資訊，說明字串樣式於. 需插入一個g於n後面，並於. 與. 各需一個替代將e換成i、將. k換成s。而其他位置的值各別表示，在其位置上字串樣式子字串與輸入文字子字串相對的edit distance。如圖2-3，. 1指出字串樣式子字串ki與. 輸入文字子字串si的edit distance為1。. 第二節. Bit-parallelism. Bit-parallelism 利用位元運算子模擬 NFA 在讀取輸入文字時，狀態啟動或關閉的結果，判斷是否有 match 發生。NFA 是指非確定有限狀態自動機，與傳統確定有限狀態自動機（Deterministic Finite Automaton, DFA）不同之處在於 NFA 允許在同一時間點時，有多個狀態啟動；而 DFA 在同一時間點時，只允許一個狀態啟動。. 10.

(27) 就上述特質而言，NFA 的特有運作方式與近似字串比對極為相似，皆可允許一些例外的情況發生。而 NFA 如同樹狀結構，需要定義資料結構與運作規則。圖 2-4 中，𝐷𝑛 的數量為字串樣式的長度﹘edit distance，即 4﹘2 ﹦2。. 圖 2-4 非確定有限狀態自動機（NFA）註：圖中為字串樣式 ATCG，輸入文字 A 之狀態. 圖 2-4 即是一非確定有限狀態自動機，以階層顯示容許之誤差值，即 edit distance；灰色代表狀態啟動，白色狀態除了在𝐷 與𝐷 之內代表狀態未啟動，其餘皆忽略其狀態。特別須注意的是，由於必須保證所有的輸入文字都會進入此 NFA，因此特別設定一個循環的機制在起點位置上；而 edit distance 由階層來顯示，故在 edit distance = 0 階層中，可忽視後兩位字元，edit distance = 1 階層則忽視前後各一位字元，edit distance = 2 則忽視前兩位字元。以圖 2-4 來看，在比對近似字串時，只須關注在 Start、 𝐷 與𝐷 兩個對角線，甚至 Start 也可忽視，因為此對角線狀態會永遠啟動。由於只需儲存𝐷 與𝐷 的變動，在實作時，以 2 進制呈現 NFA 的運作，0 表示狀態啟動，1 表示狀態未啟動。而為了避免位元運算時，造成溢位（overflows），所以兩條對角線中間必須有一個 0 當作分隔位元， 11.

(28) 圖 2-4 的狀態即是0𝐷 0𝐷 。但需特別注意𝐷 與𝐷 中數字的儲存位置，以圖 2-5 說明如下：. 圖 2-5 說明如何儲存𝐷 與𝐷. 以 0 和 1 來表示，在 bit-parallelism 之中的𝐷 為 001，𝐷 為 011，有別於一般人直覺地對應位置。而最終判別字串樣式與輸入文字是否相似的關鍵位置，在於𝐷. 011中，左起第一個位元是否為 0，其位置稱為. 最終狀態（final state）；圖 2-5 中，則其最終狀態是最右下角 3 的位置。因此，圖 2-5 的 NFA 狀態為 00010011，其中左起第一位與第五位為分隔位元，對應到每個狀態為0. 40 3. 。Bit-parallelism 之演算式，其. 中有許多位元運算子，像是 AND、OR、XOR，而 k 代表容許的 edit distance， m 則代表字串樣式的長度。 1.. ← (𝐷 ≫ (𝑘 + 2))|𝑀𝑇𝑗. 2. D ← ((𝐷 ≪ 1)|(0𝑘+ 1)𝑚. 𝑘. ). 3. & ((𝐷 ≪ (𝑘 + ))|(0𝑘+ 1)𝑚 4. & ((( + (0𝑘+ 1)𝑚 5. &(01𝑘+ )𝑚. 𝑘 )^. 𝑘. ) ≫ 1). 𝑘. 12. 01𝑘+ ).

(29) 上述各行運算式可以呈現 NFA 內部的運作，說明如下（Hyyrö, 2008）： 1.. 此行代表 NFA 水平的移動。x 記錄了𝐷. 有哪些狀態已啟動，而𝐷 需. 依據 x 來決定要啟動哪些狀態，只有𝐷. 中已啟動的狀態才可以針. 對當下讀入字元作對應的啟動。 2.. 此行代表 NFA 對角線λ的移動。向左移 1 位元，目的是檢查對角線 𝐷 裡第（𝑑 − 1）個狀態是否已啟動，若已啟動則第𝑑個狀態也將啟動。. 3.. 此行代表 NFA 垂直的移動。向左移（k+3）位元，目的是檢查𝐷 + 中第（𝑑 − 1）個狀態是否已啟動，若已啟動則啟動𝐷 中第𝑑個狀態。. 4.. 此行代表 NFA 對角線ε的移動。若𝐷 中第𝑑個狀態已啟動，則𝐷 中第（𝑑 ≤. 5.. ≤ 𝑘 + 1）個狀態也將啟動。. 此行代表遮罩，為避免溢位而設，要重新將分隔字元設為 0。在演算式中，第 2、3 行可能會造成𝐷 的第一個狀態產生錯誤，由於. 沒有前一個狀態可供檢查，故在此情況下，皆設定成前一狀態未啟動。此外，第 1 行的𝑀𝑇𝑗 代表當 match 時，NFA 應啟動的狀態，必須在比對之前預先處理。而第 3、4、5 行皆有一個由 0 與 1 組成的表示式，例如第 3 行的（0𝑘+ 1）. 𝑚 𝑘. 01𝑘+ 。此表示式之次方代表重複次數，以 m=4,. k=2 為例，此表示式為 00010111。以下將舉一簡單範例，分解每個循環 NFA 之運作。給定 k 為 2，輸入字串為 ATCGGG，字串樣式為 ATCG，則 m 為 4。首先需設定𝑀𝑇𝑗 ，如圖 2-6，𝑀𝑇𝑗 分別為𝑀𝑇1 為 01100111、𝑀𝑇2 為 01010110、𝑀𝑇3 為 00110101 及𝑀𝑇4 為 01110011。. 13.

(30) 圖 2-6 設定𝑀𝑇𝑗 遮罩. 接下來便可進入比對程序，即根據上述演算式進行運算。每當 𝐷𝑚. 𝑘. 中，由右邊算起第（𝑘 + 1）個狀態，即最終狀態 3 為啟動時，. 便回報 match，其過程如圖 2-7 所示。. 14.

(31) 圖 2-7 模擬 bit-parallelism 內部 NFA 運作流程. 第三節. 綜合分析. 本章所介紹的兩種演算法，皆有後續學者研究將其衍伸。由於傳統 dynamic programming 與 bit-parallelism 各有其優缺點，面對不同的需求時，適用的時機也不相同。Dynamic programming 可以顯而易見地發現其瓶頸在於記憶體容量過大，與字串樣式和輸入文字的長度成正比，往往在處理龐大資料量時造成記憶體容量不足；但其優點在於有充足的資訊可提供使用者做進一步的分析，例如可以了解每一個位置所採取的動作，可以協助 DNA 的對齊。反觀 bit-parallelism 記憶體容量需求小，又以位元運算為主，故速度快為此演算法的優勢；但較難得到完整的資訊，且有礙於 NFA 是以 2 進制的整數來表示，所以不能過大，此限制與字串樣式長度和 edit 15.

(32) distance 成反比關係。表 2-1 中，m 為字串樣式之長度，n 為輸入文字之長度，k 為 edit distance，以下為這兩種演算法的比較：表 2-1 Dynamic programming 與 bit-parallelism 比較表 Dynamic programming. 時間複雜度. O（𝑚𝑛） + O（𝑚 + 𝑛）. Bit-parallelism O（n）for small patterns. O（𝑚𝑘/𝜔𝑛） for longer patterns. 空間複雜度. O（𝑚𝑛）. O（（m-k）（k+2））. 字串樣式預處理. No. Yes. 回報位置. 回報起點，可得知任意位只能回報終點，無法得知起點及置之對應動作. 各位置對應動作。. 時間複雜度中，dynamic programming 在計算 edit distance 時，以需花費 O（mn），再加上 back-traverse table 也同樣需花費 O（m+n），所以整體時間需花費 O （mn）。而 bit-parallelism 因為受限於字串樣式的長度，在短字串樣式時，花費時間約 O （n），長字串樣式則需花費 O（𝑚𝑘/𝜔𝑛）。長短字串樣式花費時間的不同，是受到記憶體上的限制。以 64 位元的電腦為例，無號長整數型態可以表示 64 位數，而 NFA 大小必須不能超過 64 個狀態，即(𝑚 − 𝑘)（𝑘 + 2） ≤ 64。如果超過此限制，則需要將字串樣式以分段的方式表示。而𝜔代表電腦的位元數，mk 約等於整個 NFA 大小除以𝜔可知共需分成幾段的 NFA。. 16.

(33) 第三章. 平行處理架構. 前一章簡述了 dynamic programming 及 bit-parallelism 的運作方式與優缺點。本研究計畫透過實現於多核心 CPU 與 GPU 加速是此兩種演算法，而在進行加速之前，本章將探討平行加速的方式有哪些。以平行架構而言，dynamic programming 的可行性較 bit-parallelism 低，原因在於填表的順序與運算不獨立。反觀 bit-parallelism 在比對長度為 m 之字串樣式時，只需計算（𝑚 + 𝑘）個位置即可得知結果。因此在平行的架構下，本研究選取後者為主要加速的對象。平行處理架構分為兩種，一種為任務並行（task parallel），另一種為數據並行（data parallel）。本章分為五小節，分別為任務並行、邊界偵測的問題、數據並行 I、數據並行 II 與本章小結。說明如下：. 第一節. 任務並行. 任務並行就是將一個龐大的工作，分為許多互相獨立的任務單元，各別運算後，在統整結果之運算；如圖 3-1，原輸入字串不需分段，而每個執行緒分配不同字串樣式，各別進行比對，此種並行方式乃適用於多個字串樣式須比對的情況。而 GPU 提供使用者同一時間宣告大量的執行緒做簡短的工作，在執行緒越多的時候，效能越顯著，故 GPU 並不適用此並行模式。. 17.

(34) 圖 3-1 任務並行架構. 第二節. 邊界偵測的問題. 有別於前一節的任務並行架構，另一種並行架構即為數據並行，其中，數據並行又可以分為兩種方式，本研究稱之為數據並行 I 與數據並行 II。數據並行的平行架構，不論是數據並行 I 還是數據並行 II 都必須將整個輸入文字切斷分配給不同的執行緒同時進行比對。因此，數據並行必須特別注意區段與區段的邊界偵測問題（boundary detection problem），在每段被切割的數據中，要留意是否有字串樣式剛好被邊界所截斷，圖 3-2 即為邊界偵測問題發生的情形。. 圖 3-2 邊界偵測問題. 18.

(35) 在固定字串比對時，要防止此類問題產生，會將每條執行緒所負責的區段往後延伸 m 個位置；在比對多個不同長度的字串樣式時，則延伸其中長度最長的字串樣式之長度。然而近似字串比對比起固定字串比對又需考慮 edit distance，因此 bit-parallelism 在防止邊界偵測時，必須往後延伸（m+k）個位置，考慮是否涵蓋了容許的 edit distance 的範圍。舉例來說，給定一個長度為 m 的字串樣式，與 edit distance 為 k，則在進行近似字串比對時，若在長度（m+k）中，已 match（m+k）個字元，則可回報 match；倘若在長度（m+k）中，只 match（m-k-1）個字元，則需要繼續比對後面的位置，而在此位置上，極有可能發生邊界偵測的問題。圖 3-3 中，倘若只往後延伸（m+1）個位置，則在區段 2 的第二個位置上發生一個訊息遺漏；如果延伸（m+k）個位置，則可以避免此問題發生。. 圖 3-3 防止邊界偵測問題. 第三節. 數據並行 I. 數據並行 I 是在處理龐大資料量時，將資料量分成許多單元並分配 19.

(36) 給不同的執行緒，每條執行緒皆獨力完成自己的工作（見圖 3-4）。圖 3-4 中，將輸入字串切割分段，各條執行緒分別比對所有的字串樣式，而每條執行緒需負責的長度如 v 所示： v.. 輸入字串長度 +(總執行緒數量 – ) 總執行緒數量. + 字串樣式長度 +. 在運算式 v 中，先加上了總執行緒數量-1 後在相除，是為了避免商數為 0。相除後面加上字串樣式長度與 edit distance 是為了防止邊界問題產生。. 圖 3-4 數據並行 I. 第四節. 數據並行 II. 數據並行 II 則是給予輸入文字中每個位置一條執行緒 Xi，每個執行緒必須負責比對所有的字串樣式，而所負責的長度為（m+k）個位置。圖 3-5 即是數據並行 II 的範例。. 20.

(37) 圖 3-5 數據並行 II. 第五節. 本章小結. 採取數據並行架構，可以使一個龐大的工作，透過多執行緒分工，變成運算簡短的小工作，相較於任務平行架構，數據並行架構較適合使用於 GPU 的運作中。故本研究在 GPU 上，分別實現數據並行 I 與數據並行 II 於 bit-parallelism，分析其成效的差異。. 21.

(38) 22.

(39) 第四章. 實驗結果. 本章主要分為三個小節，分別為利用實驗環境、以 CPU 加速、以 GPU 加速及以多個 GPU 加速。說明如下：. 第一節. 實驗環境. 本實驗採用數據並行的平行架構，進行於兩台主機上，host 端與 device 端。Host 端搭配 Intel® CoreTM i7-950 CPU、Linux X86_64 作業系統、12GB DDR3 記憶體以及 ASUS P6T-SE 主機板；而 device 端則使用 NVIDIA® GeForce® GTX480 GPU 加上同樣地 CoreTM i7 系統與 NVIDIA® 驅動版本 285.05.33 和 CUDA 3.2 版，見表 4-1。而本實驗之字串樣式與輸入文字皆為 ATCG 四種字元隨機產生之文件。表 4-1 Host 與 device 主機設備規格表 host. device. Core i7 950. GTX480. Core frequency (GHz). 3.06. 1.4. Number of Streaming Multiprocessor (SM). n/a. 16. # of Cores. 4. 480. # of Threads. 8. 1,536 per SM. Single-Precision Floating Point SIMD width. 4. 32. Bandwidth GB/s. 25.6. 177.4. L1 cache, KB. 32KB per SM. 16KB per SM. On-chip shared memory KB. n/a. 48KB per SM. L2 cache, KB. 256KB per SM. 768KB for all SMs. L3 cache. 8 MB for all SM. n/a. Specification. 23.

(40) 第二節. 以 CPU 加速. 本研究於 host 主機上，利用 pthread library 產生多執行緒平行加速 dynamic programming 與 bit-parallelism。實驗結果如下：表 4-2 CPU 加速 dynamic programming 之數據 Dynamic programming 輸入文字 256Mbytes Edit distance. 執行緒數量. 1. 4. 8. 16. 32. 64. 128. 256. 512. 3. 0.1188 0.4075 0.5315 0.5056 0.5036 0.5053 0.5014 0.4473 0.3867. 5. 0.1154 0.4076 0.5079 0.5009 0.4886 0.4986 0.4967 0.4381 0.3775. 7. 0.1176 0.4075 0.5056 0.5003 0.5011 0.5031 0.5095 0.4454 0.3816. 比例. 3.4752 4.3917 4.2831 4.2447 4.2836 4.2836 3.7828 3.2569. 1. Throughput(Mbps). 4 3.5 3 2.5 Throughput (Mbps). 2. Edit distance 3. 1.5. Edit distance 5. 1. Edit distance 7. 0.5 0 1. 4. 8. 16. 32. 64. 128. 256. 512. Number of thread. 圖 4-1 以 CPU 加速 dynamic programming 效能直條圖表 4-3 CPU 加速 bit-parallelism 之數據 Bit-parallelism 輸入文字. 執行緒數量. 256Mbytes Edit distance 比例. 1. 4. 8. 16. 32. 64. 128. 256. 512. 3. 0.5806 2.3277 3.3466 3.3465 3.3452 3.3451 3.3093 3.2938 3.2993. 5. 0.5114 2.0589 3.0535 3.0526 3.0391 3.0485 3.0334 3.0243 2.9934. 7. 0.5771 2.3572 3.3932 3.3921 3.3896 3.3596 3.3584 3.3385 3.3191 1. 4.0403 5.8674 5.8661 5.8557 5.8433 5.8121 5.7855 5.7586 Throughput(Mbps) 24.

(41) 4 3.5 3 2.5 Throughput (Mbps). 2. Edit distance 3. 1.5. Edit distance 5. 1. Edit distance 7. 0.5 0 1. 4. 8. 16. 32. 64. 128. 256. 512. Number of thread. 圖 4-2 以 CPU 加速 bit-parallelism 效能直條圖. 本實驗中，字串樣式長度為 10，輸入文字長度為 256Mbytes，edit distance 以 3、5 及 7 為例。由上列資料顯示，bit-parallelism 較 dynamic programming 快上 6 倍的速度。其中， dynamic programming 在與 bit-parallelism 皆在執行緒數量為 8 時，加速達到最快。其原因在於實驗機器為四核心主機，透過超執行緒技術（hyper threading），可模擬八個邏輯處理器，正好一個邏輯處理器給與一條執行緒進行運算。然而使用 pthread library 宣告過多的執行緒，會導致效能降低。由於 dynamic programming 是透過建表來計算 edit distance，而在最後一列判斷 edit distance 的大小是否可以容許。因此無論 edit distance 為何，皆不影響其效能。另一方面，bit-parallelism 利用 NFA 的階層代表 edit distance，edit distance 越大，則階層越多，相對的儲存的值越大，位移量越多，因而效能也會下滑。另外，在實驗進行時，以 pthread library 加速 dynamic programming 容易發生記憶體不足的情況。以輸入文字為 256Mbytes，字串樣式長度為 10，edit distance 為 5 來看。倘若宣告 16 條執行緒，則每條執行緒需負責（256Mbytes＋16－1）÷16＋10＋5=32Mbytes 的區段，需建立 25.

(42) 32Mbytes×10＝320Mbytes 大小的表格，總共占據了 320Mbytes×16 的容量，這會導致程式無法執行。反觀 bit-parallelism，無論輸入文字多大，始終只需儲存 NFA 的大小，由此特性可以證明 bit-parallelism 較適合使用 GPU 來加速。. 第三節. 以 GPU 加速. 本研究以不同的執行緒數量並給予 64Mbytes 、 128Mbytes 與 256Mbytes 等不同大小的輸入文字以進行實驗。其中，根據前一節實驗結果，發現在容許誤差值為 5 的時候，效能較低，故本實驗採取此容許誤差值來實驗，而字串樣式長度同樣為 10。在前一章所介紹的數據並行平行架構中，數據並行 II 在檔案大於執行緒數量時，每條執行緒必須位移，負責下一個位置，故此方式需在宣告大量執行緒的架構下，才可得到較佳的效率。而本研究即利用 GPU 來比較數據並行 I 與數據並行 II 的效能。實驗結果如表 4-4 所示：表 4-4 GPU 加速 bit-parallelism 之數據 Bit-parallelism 輸入文字數據並行架構. 64 Mbytes. 128 Mbytes. 256 Mbytes. I. I. I. II. II. II. 64. 7.9584. 1.3289. 8.6238. 1.3738. 7.8778. 1.3915. 每個 block 中. 128. 7.1363. 1.3774. 8.8479. 1.4527. 9.2151. 1.4933. 執行緒數. 256. 5.3595. 1.3827. 7.1766. 1.4614. 8.9375. 1.5014. 512. 3.8967. 1.3855. 5.3875. 1.4637. 7.1937. 1.5053. Throughput(Mbps). 26.

(43) 10 9 8 7 6 Throughput 5 (Mbps) 4 3 2 1 0. Number of threads per block 64 Number of threads per block 128 Number of threads per block 256 Number of threads per block 512 I. II 64 MB. I. II. I. 128 MB. II. 256 MB. 圖 4-3 以 GPU 加速 bit-parallelism 效能直條圖. 由圖 4-3 可發現，數據並行 I 的效能比數據並行 II 的效能優異，其原因在於數據並行 II 中，每條執行緒負責一個輸入文字的位置，所以各執行緒都會重疊執行（m+k-1）位置。可以明顯發現，此方法與字串樣式長度與 edit distance 大小成正比關係。而數據並行 I 為了避免邊界偵測問題，只有在區段與區段間重疊執行，重疊執行的部分明顯少於數據並行 II，因此效能較高。. 第四節. 以多個 GPU 加速. 本研究進一步使用多個 GPU 透過 openMP 加速 bit-parallelism，實驗於 device 端使用兩個 NVIDIA® GeForce® GTX480 的 GPU。在實驗時，程式需宣告一個輸入文字的空間以及一個與輸入文字一樣大小的空間以儲存比對結果。而電腦的記憶體有限，不可能將所有的輸入文字讀入記憶體中，因此在龐大的輸入文字情況下，必須分段讀入輸入文字，見圖 4-4，而此次實驗以龐大輸入文字為樣本。. 27.

(44) 圖 4-4 分段讀取檔案過大的輸入文字表 4-5 Multi-GPU 加速 bit-parallelism 之數據 Bit-parallelism Single thread. 8 thread on Single GPU. on CPU. 2 GPUs. CPU. 512. 6933.11. 1858.86. 944.49. 691.15. 1024. 13901.23. 3615.57. 1886.92. 1365.02. 2048. 27971.71. 7244.54. 3779.17. 2661.74. 輸入文字 (Mbytes). Time(ms). 28.

(45) Bit-Parallelism 30000 25000 20000. Single thread on CPU. Time 15000 (ms) 10000. 8 thread on CPU Single GPU. 5000. 2 GPUs. 0 512. 1024. 2048. Input size (Mbytes). 圖 4-5 實驗結果折線圖. 本實驗隨機產生 512Mbytes、1024Mbytes 及 2048Mbytes 等龐大輸入文字分段讀入 host 端的 CPU，再由 host 端分配到 device 端的 GPU 進行比對，最後由 device 端的 GPU 複製比對結果至 host 端的 CPU 進行合併，在進行下一個輸入文字的區段。本實驗之字串樣式長度為 10，edit distance 為 5。由 4-5 的實驗結果可知，多個 GPU 處理龐大輸入文字所花費的時間較 CPU 與單個 GPU 所花費的時間少。以兩個 GPU 執行 bit-parallelism 比起單個 GPU 執行時，可減少約一半的時間。比起傳統執行於單一執行緒 CPU 版本的 bit-parallelism，在本研究 2Gbytes 的輸入文字中，以 2 個 GPU 搭配數據並行 I 的平行架構，可提高 10 倍左右的速度。. 29.

(46) 30.

(47) 第五章. 結論. 本研究透過實現 bit-parallelism 並執行於多個 GPU 上，相較於傳統執行於 CPU 上的 bit-parallelism，更能有效的加速近似字串比對。此外，本研究之 device 端為兩個 NVIDIA® GeForce® GTX480 的 GPU，由實驗結果可以發現，多個 GPU 在處理龐大資料量時，可以提供較高的效能。本研究採取 bit-parallelism 為研究對象，乃基於其速度快與節省記憶體空間，而本實驗所需結果之資訊僅判別是否有近似的字串發生。倘若所需結果之資訊為該位置上採取的動作時，並不適合採用 bit-parallelism，而較適合採取 dynamic programming 來進行比對。因此，目前仍然有許多研究學者積極想要找出速度快、省記憶體空間並擁有豐富資訊量的演算法，而透過 CPU 或 GPU 進行加速。另外，有部分學者也透過叢集運算，結合安裝在不同主機上的 CPU 或 GPU 進行運算。例如亞馬遜（Amazon）提供的雲端服務，可以給予使用者依自己需求向伺服端開啟資源或機器，將單一電腦無法於有效時間內處理完的問題，向亞馬遜要求開啟多台機器分工處理。未來研究建議透過雲端服務建置一個 GPU 叢集運算，透過連結多個不同主機上的 GPU，分別進行運算。. 31.

(48) 32.

(49) 參考文獻一、外文部分 Baeza-Yates, R., & Navarro, G. (1999). Faster approximate string matching. Algorithmica, 23(2), 127-158. Dinu, L. P., & Ionescu, R. (2011, September). A genetic approximation of closest string via rank distance. In L. Ciortuz (Chair), Artificial Intelligence III. Symposium conducted at the Symbolic and Numeric Algorithms. for. Scientific. Computing. Symbolic. and. Numeric. Algorithms for Scientific Computing, Timisoara, Romania. Dong, Y., & Qi, B. (2010, December). The technology of music retrieval by humming and its application in internet music search system. Paper presented at the 2010 IEEE International Conference on Information Theory and Information Security, Beijing, China. Hyyrö, H. (2008). Improving the bit-parallel NFA of Baeza-Yates and Navarro for approximate string matching. Information Processing Letters, 108(5), 313-319. Ivanko, E. (2006, July). Fast approximate search in strings with rearrangements. In Y. Zhong (Chair), Neural Networks. Symposium conducted at the 5th IEEE International Conference on Cognitive Informatics, Beijing, China. Kaplan, K. M., & Kaplan, J. J. (2004, October). Multiple DNA sequence approximate matching. Paper presented at the 2004 IEEE Symposium on Computational Intelligence in Bioinformatics and Computational Biology, La Jolla, CA. Kalarot, R., Morris, J., & Gimel'farb, G. (2010, November). Performance analysis of multi-resolution symmetric dynamic programming stereo on GPU. Paper presented at the 25th International Conference of Image 33.

(50) and Vision Computing New Zealand, Queenstown, New Zealand. Lyras, D. P., Sgarbas, K. N., & Fakotakis, N. D. (2007, October). Using the Levenshtein edit distance for automatic lemmatization: A case study for modern greek and english. Paper presented at the 19th IEEE International Conference on Tools with Artificial Intelligence, Paris, France. Li, H., Ni, B., Wong, M., & Leung, K. (2011, June). A fast CUDA implementation of agrep algorithm for approximate nucleotide sequence matching. Paper presented at the IEEE 9th Symposium on Application Specific Processors, San Diego, CA. Liao, H., Lin, Y., & Medioni, G. (2011, November). Aerial 3D reconstruction with line-constrained dynamic programming. Paper presented at the 13th IEEE International Conference on Computer Vision, Barcelona, Spain. Liu, W., Schmidt, B., Voss, G., Schroder, A., & Muller-Wittig, W. (2006, June). Bio-sequence database scanning on a GPU. In C. Tseng (Chair), High Performance Computational Biology. Symposium conducted at the 20th IEEE International Parallel & Distributed Processing Symposium, Rhodes Island, Greece.. Liu, T., Huang, X., Yang, L., & Zhang, P. (2009, September). Query by humming: Comparing voices to voices. Paper presented at the IEEE International Conference on Management and Service Science, Beijing, China. Liu, Y., Huang, W., Johnson, J., & Vaidya, S. (2006, May). GPU accelerated Smith-Waterman. In D. Göddeke (Chair), Database Applications Computer Graphics and Modelling. Symposium conducted at the International Conference on Computational Science, UK. Liu, Z., Lin, W., Li, N., & Lee, D. (2005, November). Detecting and filtering instant messaging spam - a global and personalized approach. Paper 34.

(51) presented at the 1st IEEE ICNP Workshop on Secure Network Protocols, Boston, MA. Myers, G. (1999). A fast bit-vector algorithm for approximate string matching based on dynamic programming. Journal of ACM, 46(3), 395-415. Moslah, O., Valles-Such, A., Guitteny, V., Couvet, S., & Philipp-Foliguet, S. (2009, May). Accelerated multi-view stereo using parallel processing capababilities of the GPUS. Paper presented at the 3DTV Conference: The True Vision - Capture, Potsdam, Germany. Munekawa, Y., Ino, F., & Hagihara, K. (2008, October). Design and implementation. of. the. Smith-Waterman. algorithm. on. the. CUDA-compatible GPU. Paper presented at the 8th IEEE International Conference on BioInformatics and BioEngineering, Athen, Greece. Qin, Z., Li, P., Zhu, Q., & Tian, C. (2010, March). SWEE: Approximately searching web service with keywords effectively and efficiently. Paper presented at the 2nd IEEE International Conference on Advanced Computer Control, Shenyang, China. Si, J., Yang, L., Lu, C., Sun, J., & Mei, S. (2009, June). Approximate dynamic programming for continuous state and control problems. Paper presented at the 17th Mediterranean Conference on Control and Automation, Thessaloniki, Greece. Soleh, M. Y., & Purwarianti, A. (2011, July). A non word error spell checker for Indonesian using morphologically analyzer and HMM. Paper presented at the 2011 IEEE International Conference on Electrical Engineering and Informatics, Bandung, Indonesia. Smith, T. F., & Waterman, M. S. (1981). Identification of common molecular subsequences. Journal of Molecular Biology, 147(1), 195-197. Song, T., Xue, Y., & Wang, D. (2006, October). An algorithm of large-scale 35.

(52) approximate multiple string matching for network security. In H. Chen (Chair), Security Protocols and Watermarks. Symposium conducted at the Communications and Networking in China, Beijing, China. Shi, W., & Xie, M. (2011, June). Spam filtering cloud platform based on sharing fingerprints. Paper presented at the 2011 IEEE International Conference on Computer Science and Service System, Nanjing, China. Wu, O., Zuo, H., Hu, W., Zhu, M., & Li, S. (2008, December). Recognizing and filtering web images based on people's existence. Paper presented at the IEEE/WIC/ACM International Conference on Web Intelligence and Intelligent Agent Technology, Sydney, Australia. Wu, S., & Manber, U. (1992). Fast text searching: allowing errors. COMMUNICATION OF THE ACM, 35(10), 83-91.. 36.

(53) 附. 錄. 37.

(54) 38.

(55) 附錄一、Dynamic programming（CPU 版本） #include <stdio.h> #include <stdlib.h> #include <string.h> #include <time.h> #include <math.h> #include <pthread.h> #define PATTERN_LEN 99 #define NUM_THREADS 1 #define iMIN(x,y) ((x)>(y))? (y):(x) pthread_t tid[NUM_THREADS]; pthread_attr_t attr; void *parallel_bitap_NFA(void *param); int value[NUM_THREADS]={0}; void *fill_dp(void *param); static int pattern_num = 0; unsigned long int pattern_code = 0; int k; int dp_width; int dp_depth; int pattern_len, input_len; char* input; char pattern[100]; unsigned long int *dp_result; unsigned long int Temp1 = 0; static int dp_count = 0; int main(int argc, char* argv[]){ FILE* fr_input; FILE* fr_pattern; FILE* fw_dp_match; int i,j,y,policy; int bit_length = 64; struct timeval t_start, t_end; float elapsedTime; if(argc!=4){ printf("using command:%s input pattern k\n",argv[0]); return 0; } if((fr_input=fopen(argv[1],"rb"))==NULL) { printf("Input FILE: %s cannot be opened\n",argv[1]); exit(1); } /* obtain file size: */ 39.

(56) fseek (fr_input , 0, SEEK_END); input_len = ftell (fr_input); printf("original input_len is %d\n", input_len); rewind(fr_input); /* allocate memory to contain the whole file: */ input = (char *) malloc (sizeof(char)*input_len); if (input == NULL) { printf("Allocate input memory error"); exit (1); } /* copy the file into the buffer: */ if(input_len != fread (input, 1, input_len, fr_input)){ printf("error occurs when reading input file to memory"); exit(1); } if(input[input_len-1]=='\n'){ input[input_len-1]='\0'; input_len=strlen(input); printf("exact input_len is %d\n", input_len); } fclose(fr_input); //malloc result array dp_result and dp2m_result dp_result = (unsigned long int*)malloc(sizeof(unsigned long int)*input_len); for(j=0;j<input_len+1;j++){ dp_result[j] = 0; } //reading pattern file if((fr_pattern=fopen(argv[2],"rb"))==NULL) { printf("Pattern FILE: %s cannot be opened\n",argv[2]); exit(1); } k=atoi(argv[3]);//reading edit distance /**********************Pthread***********************************/ /* get the default attributes */ pthread_attr_init(&attr); /* get the current scheduling policy */ if (pthread_attr_getschedpolicy(&attr,&policy) != 0) fprintf(stderr, "Unable to get policy.\n"); else { 40.

(57) if (policy == SCHED_OTHER) printf("SCHED_OTHER\n"); else if (policy == SCHED_RR) printf("SCHED_OTHER\n"); else if (policy == SCHED_FIFO) printf("SCHED_FIFO\n"); } /* set the scheduling policy - FIFO, RT, or OTHER */ if (pthread_attr_setschedpolicy(&attr, SCHED_OTHER) != 0) printf("unable to set scheduling policy to SCHED_OTHER \n"); printf("*******************DP Method********************\n"); for(y=0;y<NUM_THREADS;y++) value[y]=y; // start time gettimeofday(&t_start, NULL); //dp function while(fscanf(fr_pattern, "%s", pattern) == 1 ){ pattern_code = pow(2, pattern_num); pattern_len=strlen(pattern); if(pattern[pattern_len-1]=='\n'){ pattern[pattern_len-1]='\0'; pattern_len=strlen(pattern); printf("pattern_len is %d\n", pattern_len); } dp_width=input_len+1; dp_depth=pattern_len+1; for (y = 0; y < NUM_THREADS; y++){ pthread_create(&tid[y],&attr,fill_dp,(void *)&value[y]); } for (y = 0; y < NUM_THREADS; y++) pthread_join(tid[y], NULL); pattern_num++; } // stop time gettimeofday(&t_end, NULL); /* compute and print the elapsed time in millisec */ elapsedTime = (t_end.tv_sec - t_start.tv_sec) * 1000.0; elapsedTime += (t_end.tv_usec - t_start.tv_usec) / 1000.0; printf("The input size is %d\n", input_len); printf("The elapsed time is %f ms\n", elapsedTime); printf("The throughput is %f Mps\n",(float)(input_len*8)/(elapsedTime*1000000) ); fclose(fr_pattern); // Output matching results 41.

(58) fw_dp_match = fopen("dp_approximate_match_result.txt", "w"); if (fw_dp_match == NULL) { perror("Open output file failed.\n"); exit(1); } for(i=0;i<input_len+1;i++){ if(dp_result[i] != 0){ for (j = pattern_num; j >= 0; j--) { Temp1 = (dp_result[i] >> j) & 1; if (Temp1 == 1) { (void)fprintf(fw_dp_match, "At position %4d, match pattern %d\n", i+1, j+1); dp_count++; } } } } fclose(fr_pattern); fclose(fw_dp_match); free(dp_result); return 0; } void *fill_dp(void *param){ int segmentSize=(input_len+NUM_THREADS-1)/NUM_THREADS; int remainder=pattern_len+k; int tid=(*(int*)param); int pos=segmentSize*tid ; int end=pos+segmentSize+remainder; end=iMIN(end,input_len); int *dp; int distance=end-pos+1; dp = (int *)malloc(distance*dp_depth*sizeof(int)); initial_dp(dp, distance, dp_depth); int i; int j; for(i=1;i<dp_depth; i++){ for(j=pos+1;j<end+1; j++){ if(input[j-1]==pattern[i-1]) dp[i*distance+(j-pos)] = minThree(dp[(i-1)*distance+(j-pos)]+1,dp[i*distance+(j-pos)-1]+1,dp[(i-1)*distance +(j-pos)-1]); else{ dp[i*distance+(j-pos)] = minThree(dp[(i-1)*distance+(j-pos)]+1,dp[i*distance+(j-pos)-1]+1,dp[(i-1)*distance 42.

(59) +(j-pos)-1]+1); } } } for(i=pos+1; i<end+1;i++){// while the variable end add 1 for the reason of initial entry if(dp[(dp_depth-1)*distance+(i-pos)]<=k) dp_result[i-1] |= pattern_code; } free(dp); } void initial_dp(int* dp_ptr, int width, int depth){ int i; for(i=0;i<width;i++) dp_ptr[i]=0; for(i=1;i<depth;i++) dp_ptr[i*width]=i; }. 43.

(60) 附錄二、Bit-parallelism（CPU 版本） #include <stdio.h> #include <stdlib.h> #include <string.h> #include <limits.h> #include <math.h> #include <sys/time.h> #include <pthread.h> #define PATTERN_LEN 99 #define iMIN(x,y) ((x)>(y))? (y):(x) #define NUM_THREADS 1 pthread_t tid[NUM_THREADS]; pthread_attr_t attr; void *parallel_bitap_NFA(void *param); int value[NUM_THREADS]={0}; char *text; char pattern[100]; int k, bit_length; unsigned long int text_size = 0; static int m; unsigned long int *match_result; unsigned long int pattern_code = 0; static int pattern_num = 0; struct timeval t_start, t_end; float elapsedTime; int main (int argc, const char *argv[]) { FILE *fpin; FILE *fr; FILE *fpout; int y, r, policy;; static int count = 0; unsigned long int Temp = 0; if (argc != 4) { printf("command: %s inputFile pattern distance\n", argv[0]); exit(1); } if ((fpin = fopen(argv[1],"rb")) == NULL) { printf("Input FILE: %s cannot be opened\n", argv[1]); exit(2); } /* obtain file size: */ fseek(fpin , 0 , SEEK_END); text_size = ftell(fpin); 44.

(61) rewind(fpin); /* allocate memory to contain the whole file: */ text = (char*) malloc (sizeof(char) * text_size); if (text == NULL) { perror("Allocate input memory error"); exit (3); } /* copy the file into the buffer: */ y = fread (text, 1, text_size, fpin); if (y == 0){ printf("Copy fail!\n"); exit(4); } fclose(fpin); /* allocate memory for output */ match_result = (unsigned long int *) malloc (sizeof(unsigned long int) * text_size); if (match_result == NULL) { perror("Allocate output memory error"); exit (5); } for (y = 0; y < text_size; y++) { match_result[y] = 0; } //asign pattrtn if((fr = fopen(argv[2],"rb")) == NULL) { printf("Pattern FILE: %s cannot be opened\n", argv[2]); exit (6); } // copy distance from command line k = atoi(argv[3]); /**********************Pthread***********************************/ /* get the default attributes */ pthread_attr_init(&attr); /* get the current scheduling policy */ if (pthread_attr_getschedpolicy(&attr,&policy) != 0) fprintf(stderr, "Unable to get policy.\n"); else { if (policy == SCHED_OTHER) printf("SCHED_OTHER\n"); else if (policy == SCHED_RR) printf("SCHED_OTHER\n"); else if (policy == SCHED_FIFO) 45.

(62) printf("SCHED_FIFO\n"); } /* set the scheduling policy - FIFO, RT, or OTHER */ if (pthread_attr_setschedpolicy(&attr, SCHED_OTHER) != 0) printf("unable to set scheduling policy to SCHED_OTHER \n"); for(y=0;y<NUM_THREADS;y++) value[y]=y; // start time gettimeofday(&t_start, NULL); /* call the algorithm */ while (fscanf(fr, "%s", pattern) == 1 ) { puts(pattern); m = strlen(pattern); for (y = 0; y < NUM_THREADS; y++){ pthread_create(&tid[y],&attr,parallel_bitap_NFA,(void *)&value[y]); } /* Now join on each thread */ for (y = 0; y < NUM_THREADS; y++) pthread_join(tid[y], NULL); pattern_num++; } // stop time gettimeofday(&t_end, NULL); /* compute and print the elapsed time in millisec */ elapsedTime = (t_end.tv_sec - t_start.tv_sec) * 1000.0; elapsedTime += (t_end.tv_usec - t_start.tv_usec) / 1000.0; printf("The input size is %ld bytes\n", text_size ); printf("The elapsed time is %lf ms\n", elapsedTime); printf("The throughput is %lf Mbps\n",((double)(text_size*8))/(elapsedTime*1000000)); // Output results fpout = fopen("Bitap_Pthread_match_result.txt", "w"); if (fpout == NULL) { perror("Open output file failed.\n"); exit(7); } // Output match result to file for (y = 0; y < text_size; y++) { if (match_result[y] != 0) { for (r = 63; r >= 0; r--) { Temp = (match_result[y] >> r) & 1; if (Temp == 1) { fprintf(fpout, "At position %4d, match pattern %d\n", y+1, r+1); 46.

(63) count++; } } } } fclose(fpout); fclose(fr); free(text); free(match_result); return 0; } void *parallel_bitap_NFA(void *param){ int segmentSize = (text_size + NUM_THREADS-1) / NUM_THREADS ; int remainder = m + k; int tid=(*(int*)param); int pos = segmentSize * tid ; int end = pos + segmentSize + remainder; end = iMIN(end, text_size) ; int i, j, z; int *Mtable; unsigned long D ; unsigned long X = 0; unsigned long *M; unsigned long mask1 = 0; unsigned long mask2 = 0; unsigned long mask3 = 0;. bit_length = (k+2) * (m-k); if (bit_length > 64) { perror("pattern is larger than 64 bits!"); exit(8); } if (pattern[0] == '\0') { printf("No pattern input!"); exit(9); } pattern_code = pow(2, pattern_num); Mtable = (int*)malloc(sizeof(int)* m * m); /*Construct table in order to set up mask*/ for (i = 0; i < m; i++) { for (j = 0; j < m; j++) { if (pattern[i] == pattern[j]) { Mtable[i*m+j] = 0; //equal 0 }else { Mtable[i*m+j] = 1; //nonequal 1 47.

(64) } } } /* Initialize the bit array */ M = malloc(sizeof(unsigned long int) * m); for (i = 0; i < m; i++) M[i] = 0; /*Set up the mask of each element of the pattern*/ for (i = 0; i < m; i++) { z = 0; for (j = 0; j < bit_length; j++) { if ((j+1)%(k+2) == 0) { z++; //shift D continue; }else { M[i] += (Mtable[i*m+m-k-1+(j%(k+2))-z] * pow(2, j)); } } } /*Set up the mask using throughout the algorithm*/ /**mask 1**/ for (i = 0; i < bit_length; i++) { if (i%(k+2) == 0) { mask1 += pow(2, i); } } /**mask 2**/ for (i = 0; i < bit_length; i++) { if (i < k+1) { mask2 += pow(2, i); }else { if (i%(k+2) == 0) { mask2 += pow(2, i); } } } /**mask 3**/ for (i = 0; i < bit_length; i++) { if ((i+1)%(k+2) != 0) mask3 += pow(2, i); } /*Strat matching*/ D = mask3; 48.

(65) while (pos < end) { j = 0; while (j < m) { if (text[pos] == pattern[j]) { X = ((D >> k) >> 2) | M[j] ;/*Line 1*/ break; } j++; } if (j == m) X = ((D >> k) >> 2) | mask3 ; D = ((D << 1) | mask1) & (((D << k) << 3) | mask2) & (((X + mask1) ^ X) >> 1) & mask3 ;/*Line 2 - Line 5*/ if(((D >> k) & 1) == 0 ) { match_result[pos] |= pattern_code; } pos++; } free(M); free(Mtable); }. 49.

(66) 附錄三、Bit-parallelism（GPU 版本） #include <stdio.h> #include <stdlib.h> #include <string.h> #include <limits.h> #include <math.h> #include <sys/time.h> #include <cuda_runtime.h> #define MAX_MB 512 * 1024 * 1024 #define iMIN(x,y) ((x)>(y))? (y):(x) void parallel_bitap_NFA(char *match_result, unsigned long int text_size, int pattern_number ); __global__ void kernel( char *text, char *pattern, int segmentSize, int Overlap, unsigned long int *M, char *match_result, int m, int text_size, int k, unsigned long int mask1, unsigned long int mask2, unsigned long int mask3, char pattern_code ); int k, m, bit_length ; char *input ; char pattern[100] ; char *m_result ; cudaError_t status; unsigned long int input_len; unsigned long int seg_text_length = 0; unsigned long int seg_text_start = 0; unsigned long int seg_text_end = 0; double TotalsetupInputTime = 0; double TotalgetResultTime = 0; double GPUTime = 0; int main (int argc, const char *argv[]) { FILE *fr_input; FILE *fr_pattern; FILE *fpout; int i; int count = 0; int pattern_num; size_t result; if (argc != 4) { printf("command: %s inputFile pattern distance\n", argv[0]); exit(1); } if ((fr_input = fopen(argv[1],"rb")) == NULL) { printf("Input FILE: %s cannot be opened\n", argv[1]); exit(1); } 50.

(67) //asign pattrtn if((fr_pattern = fopen(argv[2],"rb")) == NULL) { printf("Pattern FILE: %s cannot be opened\n", argv[2]); exit (1); } // copy distance from command line k = atoi(argv[3]); // Output results fpout = fopen("Bitap_LargeText_GPU_match_result.txt", "w"); if (fpout == NULL) { perror("Open output file failed.\n"); exit(1); } /********For GPU device********/ int deviceID = 0; if ( cudaSetDevice(deviceID) != cudaSuccess ){ fprintf(stderr, "Set CUDA device %d error\n", deviceID); exit(1); } cudaDeviceProp deviceProp; cudaGetDeviceProperties(&deviceProp, deviceID); printf("Using Device %d: \"%s\"\n", deviceID, deviceProp.name); /******************************/ printf("*******************Bitap Method on GPU********************\n"); /* obtain file size: */ fseek(fr_input , 0 , SEEK_END); input_len = ftell(fr_input); printf("original input_len is %ld\n", input_len); rewind(fr_input); if ( input_len <= MAX_MB) { seg_text_length = input_len; }else{ seg_text_length = MAX_MB; } seg_text_end = seg_text_length - 1 ; while( seg_text_end < input_len) { pattern_num = 1; seg_text_length = seg_text_end - seg_text_start + 1; input = (char *) malloc (sizeof(char)*seg_text_length); if (input == NULL){ if (seg_text_end != 0){ break; }else{ printf("Allocate input memory error\n"); 51.

(68) exit (1); } } /* copy the file into the buffer: */ fseek (fr_input , seg_text_start, SEEK_SET); result = fread (input, 1, seg_text_length, fr_input); if( result != seg_text_length ){ printf("error occurs when reading input file to memory\n"); exit(1); } //malloc result array dp_result and dp2m_result m_result = (char*) malloc (sizeof(char) * seg_text_length); for(i=0;i<seg_text_length;i++){ m_result[i] = 0; } //call the algorithm while (fscanf(fr_pattern, "%s", pattern) == 1 ) { puts(pattern); m = strlen(pattern); bit_length = (k+2) * (m-k); parallel_bitap_NFA(m_result, seg_text_length, pattern_num); pattern_num++; } // Output match result to file for (i = 0; i < seg_text_length; i++) { if (seg_text_start > 0 && i < (m + k)) continue; if (m_result[i] > 0) { fprintf(fpout, "At position %4ld, match pattern %d\n", i+seg_text_start+1, pattern_num-1); count++; } } seg_text_start = seg_text_start + MAX_MB - m - k; seg_text_end = seg_text_start + MAX_MB - 1; seg_text_end = iMIN(seg_text_end, input_len-1); rewind(fr_input); rewind(fr_pattern); } printf("\nThe gpu elapsed time is %.2f ms\n", GPUTime); printf("The input length is %.2f Mbytes\n", (float)input_len/1024/1024); printf("The throughput is %.2f Mbps\n",(float)(input_len*8)/(GPUTime*1000000/8) ); printf("\nsetupInputTime = %.2f ms\n", TotalsetupInputTime); 52.

(69) printf("getResultTime = %.2f ms\n", TotalgetResultTime); double totalTime = GPUTime + TotalsetupInputTime + TotalgetResultTime ; printf("total time including PCIe = %.2f ms\n", totalTime); printf("The throughput(PCIe) is %lf Mbps\n",(float)(input_len*8)/(totalTime*1000000) ); printf("The bandwidth of CPU to GPU is %.2f Gbps\n",(float)(input_len)/(TotalsetupInputTime*1000000/8) ); printf("The bandwidth of GPU to CPU is %.2lf Gbps\n",(double)(input_len)/(TotalgetResultTime*1000000/32) ); fclose(fpout); fclose(fr_pattern); free(input); free(m_result); return 0; } void parallel_bitap_NFA(char *d_match_result, unsigned long int text_size, int pattern_num ){ struct timeval t_start, t_end; float times; unsigned long int mask1 = 0; unsigned long int mask2 = 0; unsigned long int mask3 = 0; unsigned long int *M; if (bit_length > 64) { perror("pattern is larger than 64 bits!"); exit(1); } int i, j, z; int *Mtable; if (pattern[0] == '\0') { printf("No pattern input!"); exit(1); } /********For GPU device********/ char *dev_text; char *dev_match_result; char *dev_pattern; unsigned long int *dev_M; //Allocate memory in GPU status = cudaMalloc((void **) &dev_text, sizeof(char) * text_size); if ( cudaSuccess != status ){ fprintf(stderr, "cudaMallocHost dev_text error: %s\n", cudaGetErrorString(status)); exit(1) ; 53.

(70) } status = cudaMalloc((void **) &dev_match_result, sizeof(char) * text_size); if ( cudaSuccess != status ){ fprintf(stderr, "cudaMallocHost dev_match_result error: %s\n", cudaGetErrorString(status)); exit(1) ; } cudaMemset(dev_match_result, 0, text_size * sizeof(char)); status = cudaMalloc((void **) &dev_pattern, sizeof(char) * m); if ( cudaSuccess != status ){ fprintf(stderr, "cudaMallocHost dev_pattern error: %s\n", cudaGetErrorString(status)); exit(1) ; } status = cudaMalloc((void **) &dev_M, sizeof(unsigned long int) * m); if ( cudaSuccess != status ){ fprintf(stderr, "cudaMallocHost dev_M error: %s\n", cudaGetErrorString(status)); exit(1) ; } /******************************/ Mtable = (int*)malloc(sizeof(int)* m * m); M = (unsigned long int*)malloc(sizeof(unsigned long int) * m); for (i = 0; i < m; i++) M[i] = 0; gettimeofday(&t_start, NULL); /*Construct table in order to set up mask*/ for (i = 0; i < m; i++) { for (j = 0; j < m; j++) { if (pattern[i] == pattern[j]) { Mtable[i*m+j] = 0; //equal 0 }else { Mtable[i*m+j] = 1; //nonequal 1 } } } /*Set up the mask of each element of the pattern*/ for (i = 0; i < m; i++) { z = 0; for (j = 0; j < bit_length; j++) { if ((j+1)%(k+2) == 0) { z++; //shift D continue; }else { 54.

(71) M[i] += (Mtable[i*m+m-k-1+(j%(k+2))-z] * pow(2, j)); } } } /*Set up the mask using throughout the algorithm*/ /**mask 1**/ for (i = 0; i < bit_length; i++) { if (i%(k+2) == 0) { mask1 += pow(2, i); } } /**mask 2**/ for (i = 0; i < bit_length; i++) { if (i < k+1) { mask2 += pow(2, i); }else { if (i%(k+2) == 0) { mask2 += pow(2, i); } } } /**mask 3**/ for (i = 0; i < bit_length; i++) { if ((i+1)%(k+2) != 0) mask3 += pow(2, i); } /*************For GPU device*************/ //Copt text from Host to Device status = cudaMemcpy(dev_text, input, sizeof(char) * text_size, cudaMemcpyHostToDevice); if ( cudaSuccess != status ){ fprintf(stderr, "Error: %s\n", cudaGetErrorString(status)); exit(1) ; } status = cudaMemcpy(dev_pattern, pattern, sizeof(char) * m, cudaMemcpyHostToDevice); if ( cudaSuccess != status ){ fprintf(stderr, "Error: %s\n", cudaGetErrorString(status)); exit(1) ; } status = cudaMemcpy(dev_M, M, sizeof(unsigned long int) * m, cudaMemcpyHostToDevice); if ( cudaSuccess != status ){ fprintf(stderr, "Error: %s\n", cudaGetErrorString(status)); exit(1) ; } 55.

(72) gettimeofday(&t_end, NULL); // compute and print the elapsed time in millisec double setupInputTime = (t_end.tv_sec - t_start.tv_sec) * 1000.0; setupInputTime += (t_end.tv_usec - t_start.tv_usec) / 1000.0; TotalsetupInputTime += setupInputTime; /*****************DONE*******************/ // set 256 threads per block, set grid size automatically int dimBlock = 128 ; int dimGrid ; dimGrid = (text_size + dimBlock - 1) / dimBlock ; if (dimGrid > 65535) { dimGrid = 65535 ; } int segmentSize = (text_size + dimBlock * dimGrid - 1) / (dimBlock * dimGrid); int segmentOverlap = segmentSize + m + k ; // record time setting cudaEvent_t start, stop; cudaEventCreate(&start); cudaEventCreate(&stop); cudaEventRecord(start, 0); /*****************Call the Kernel function*****************/ kernel <<< dimGrid, dimBlock >>>( dev_text, dev_pattern, segmentSize, segmentOverlap, dev_M, dev_match_result, m, text_size, k, mask1, mask2, mask3, pattern_num ); // record time setting cudaEventRecord(stop, 0); cudaEventSynchronize(stop); cudaEventElapsedTime(&times, start, stop); cudaEventDestroy(start); cudaEventDestroy(stop); GPUTime += times ; /********Copy result form Device to Host********/ gettimeofday(&t_start, NULL); // copy the result data from device to host status = cudaMemcpy(d_match_result, dev_match_result, sizeof(char) * text_size, cudaMemcpyDeviceToHost); if ( cudaSuccess != status ){ fprintf(stderr, "Error: %s\n", cudaGetErrorString(status)); exit(1) ; } gettimeofday(&t_end, NULL); // compute and print the elapsed time in millisec double getResultTime = (t_end.tv_sec - t_start.tv_sec) * 1000.0; getResultTime += (t_end.tv_usec - t_start.tv_usec) / 1000.0; 56.

(73) TotalgetResultTime += getResultTime; /******************************/ free(M) ; cudaFree(dev_M) ; cudaFree(dev_pattern) ; cudaFree(dev_text); cudaFree(dev_match_result); } __global__ void kernel(char *text, char *pattern, int segmentSize, int Overlap, unsigned long int *M, char *dev_match_result, int m, int input_size, int k, unsigned long int mask1, unsigned long int mask2, unsigned long int mask3, char pattern_code ){ int j ; int tid = threadIdx.x + blockIdx.x * blockDim.x ; int pos = tid * segmentSize ; int end = pos + Overlap; end = iMIN(end, input_size) ; unsigned long int X = 0 ; unsigned long int D = mask3 ; while ( pos < end ){ j=0; while (j < m) { if (text[pos] == pattern[j]) { X = ((D >> k) >> 2) | M[j] ; //Line 1 break ; } j++ ; } if (j == m) X = ((D >> k) >> 2) | mask3 ; D = ((D << 1) | mask1) & (((D << k) << 3) | mask2) & (((X + mask1) ^ X) >> 1) & mask3 ;//Line 2 - Line 5 if(((D >> k) & 1) == 0 ){ dev_match_result[pos] |= pattern_code; } pos++ ; } }. 57.

(74) 附錄四、Bit-parallelism（Multi-GPU 版本） #include <omp.h> #include <stdio.h> #include <stdlib.h> #include <string.h> #include <limits.h> #include <math.h> #include <sys/time.h> #include <cuda_runtime.h> #define MAX_MB 512 * 1024 * 1024 #define BASIC_UNIT 128 * 1024 * 1024 #define MAX_THREADS 8 #define iMIN(x,y) ((x)>(y))? (y):(x) void parallel_bitap_NFA(char *d_input, char *match_result, unsigned long int text_size, int pattern_number , int num_GPU); __global__ void kernel( char *text, char *pattern, int segmentSize, int Overlap, unsigned long int *M, int shiftarray, char *match_result, int m, int text_size, int k, unsigned long int mask1, unsigned long int mask2, unsigned long int mask3, char pattern_code ); int k, m, bit_length ; char *input ; char pattern[100] ; char *m_result ; cudaError_t status; unsigned long int input_len; unsigned long int seg_text_length = 0; unsigned long int seg_text_start = 0; unsigned long int seg_text_end = 0; double TotalsetupInputTime = 0; double TotalgetResultTime = 0; double GPUTime = 0; int main (int argc, const char *argv[]) { FILE *fr_input; FILE *fr_pattern; FILE *fpout; int i; int count = 0; int pattern_num; size_t result; if (argc != 4) { printf("command: %s inputFile pattern distance\n", argv[0]); exit(1); } if ((fr_input = fopen(argv[1],"rb")) == NULL) { 58.

(75) printf("Input FILE: %s cannot be opened\n", argv[1]); exit(1); } //asign pattrtn if((fr_pattern = fopen(argv[2],"rb")) == NULL) { printf("Pattern FILE: %s cannot be opened\n", argv[2]); exit (1); } // copy distance from command line k = atoi(argv[3]); // Output results fpout = fopen("Bitap_omp_GPU_match_result.txt", "w"); if (fpout == NULL) { perror("Open output file failed.\n"); exit(1); } printf("*******************Bitap Method on GPU using OpenMP********************\n"); /* obtain file size: */ fseek(fr_input , 0 , SEEK_END); input_len = ftell(fr_input); printf("original input_len is %ld\n", input_len); rewind(fr_input); if ( input_len <= MAX_MB) { seg_text_length = input_len; }else{ seg_text_length = MAX_MB; } seg_text_end = seg_text_length - 1 ; int max_gpus = 0; int num_gpus = 0; int basic_unit = BASIC_UNIT; cudaGetDeviceCount(&max_gpus); if( 1 > max_gpus ) { printf("Error: no CUDA capable devices were detected\n"); exit(1); } int step = 1; while( seg_text_end < input_len) { pattern_num = 1; seg_text_length = seg_text_end - seg_text_start + 1; input = (char *) malloc (sizeof(char)*seg_text_length); if (input == NULL){ if (seg_text_end != 0){ 59.