實現高性能糾刪碼引擎 | 糾刪碼技術詳解（下）

作者介紹： 

徐祥曦，七牛云工程師，獨立開發了多套高性能糾刪碼/再生碼編碼引擎。
柳青，華中科技大學博士，研究方向為基于糾刪碼的分布式存儲系統。

前言：

在上篇《如何選擇糾刪碼編碼引擎》中，我們簡單了解了 Reed-Solomon Codes（RS 碼）的編／解碼過程，以及編碼引擎的評判標準。但并沒有就具體實現進行展開，本篇作為《糾刪碼技術詳解》的下篇，我們將主要探討工程實現的問題。

這里先簡單提煉一下實現高性能糾刪碼引擎的要點：首先，根據編碼理論將矩陣以及有限域的運算工程化，接下來主要通過 SIMD 指令集以及緩存優化工作來進行加速運算。也就是說，我們可以將 RS 的工程實現劃分成兩個基本步驟：

1. 將數學理論工程化

2. 進一步的工程優化

這需要相關研發工程師對以下內容有所掌握：

1. 有限域的基本概念，包括有限域的生成與運算

2. 矩陣的性質以及乘法規則

3. 計算機體系結構中關于 CPU 指令以及緩存的理論

接下來，我們將根據這兩個步驟并結合相關基礎知識展開實現過程的闡述。

一 、理論工程化

以 RS 碼為例，糾刪碼實現于具體的存儲系統可以分為幾個部分：編碼、解碼和修復過程中的計算都是在有限域上進行的；編碼過程即是計算生成矩陣（范德蒙德或柯西矩陣）和所有數據的乘積；解碼則是計算解碼矩陣（生成矩陣中某些行向量組成的方陣的逆矩陣）和重建數據的乘積。

1.1 有限域運算

有限域是糾刪碼中運算的基礎域，所有的編解碼和重建運算都是基某個有限域的。不止是糾刪碼，一般的編碼方法都在有限域上進行，比如常見的AES加密中也有有限域運算。使用有限域的一個重要原因是計算機并不能精確執行無限域的運算，比如有理數域和虛數域。

此外，在有限域上運算另一個重要的好處是運算后的結果大小在一定范圍內，這是因為有限域的封閉性決定的，這也為程序設計提供了便利。比如在 RS 中，我們通常使用 GF(2^8)，即 0～255 這一有限域，這是因為其長度剛好為1字節，便于我們對數據進行存儲和計算。

在確定了有限域的大小之后，通過有限域上的生成多項式可以找到該域上的生成元^[1]，進而通過生成元的冪次遍歷有限域上的元素，利用這一性質我們可以生成相應的指數表。通過指數表我們可以求出對數表，再利用指數表與對數表最終生成乘法表。關于本原多項式的生成以及相關運算表的計算可以參考我在開源庫中的數學工具。[2]

有了乘法表，我們就可以在運算過程中直接查表獲得結果，而不用進行復雜的多項式運算了。同時也不難發現，查表優化將會成為接下來工作的重點與難點。

1.2 選擇生成矩陣

生成矩陣(GM, generator matrix) 定義了如何將原始數據塊編碼為冗余數據塊，RS 碼的生成矩陣是一個 n 行 k 列矩陣，將 k 塊原始數據塊編碼為 n 塊冗余數據塊。如果對應的編碼是系統碼(比如 RAID)，編碼后包含了原始數據，則生成矩陣中包含一個 k×k 大小的單位矩陣和(nk)×k 的冗余矩陣， 單位矩陣對應的是原始數據塊，冗余矩陣對應的是冗余數據塊。非系統碼沒有單位矩陣，整個生成矩陣都是冗余矩陣，因此編碼后只有冗余數據塊。通常我們會使用系統碼以提高數據提取時的效率，那么接下來我們需要找到合適的冗余矩陣。

在解碼過程中我們要對矩陣求逆，因此所采用的矩陣必須滿足子矩陣可逆的性質。目前業界應用最多的兩種矩陣是 Vandermonde matrix (范德蒙矩陣）和Cauchy matrix（柯西矩陣）。其中范德蒙矩陣歷史最為悠久，但需要注意的是我們并不能直接使用范德蒙矩陣

作為生成矩陣，而需要通過高斯消元后才能使用，這是因為在編碼參數（k+m)比較大時會存在矩陣不可逆的風險。

柯西矩陣運算簡單，只不過需要計算乘法逆元，我們可以提前計算好乘法逆元表以供生成編碼矩陣時使用。創建以柯西矩陣為生成矩陣的編碼矩陣的偽代碼如下圖所示：

*// m 為編碼矩陣*

*// rows為行數，cols為列數*

*// *k×k 的單位矩陣

**for **j := 0; j < cols; j++ {

   m[j][j] = byte(1)

}

*//* mxk 的柯西矩陣

**for **i := cols; i < rows; i++ {

   **for **j := 0; j < cols; j++ {

      d := i ^ j

      a := inverseTable[d]    *// 查乘法逆元表*

      m[i][j] = byte(a)

   }

}

 1.3 矩陣求逆運算

有限域上的求逆方法和我們學習的線性代數中求逆方法相同，常見的是高斯消元法，算法復雜度是 O(n^3)。過程如下：

1. 在待求逆的矩陣右邊拼接一個單位矩陣

2. 進行高斯消元運算

3. 取得到的矩陣左邊非單位矩陣的部分作為求逆的結果，如果不可逆則報錯。

我們在實際的測試環境中發現，矩陣求逆的開銷還是比較大的(大約 6000 ns/op)。考慮到在實際系統中，單盤數據重建往往需要幾個小時或者更長（磁盤I/O 占據絕大部分時間），求逆計算時間可以忽略不計。

二、進一步的工程優化

2.1 利用 SIMD 加速有限域運算

從上一篇文章可知，有限域上的乘法是通過查表得到的，每個字節和生成矩陣中元素的乘法結果通過查表得到，圖1 給出了按字節對原始數據進行編碼的過程（生成多項式為x^8 + x^4 + x^3 + x^2 + 1）。
對于任意 1 字節來說，在 GF(2^8) 內有256種可能的值，所以沒有元素對應的乘法表大小為 256 字節。每次查表可以進行一個字節數據的乘法運算，效率很低。

圖 1, 按字節對原始數據進行編碼

目前主流的支持 SIMD 相關指令的寄存器有 128bit（XMM 指令)、256bit (YMM 指令)這兩種容量，這意味著對于64位的機器來說，分別提供了2到4倍的處理能力，我們可以考慮采用 SIMD 指令并行地為更多數據進行乘法運算。

但每個元素的乘法表的大小為 256 Byte ,這大大超出了寄存器容納能力。為了達到利用并行查表的目的，我們采用分治的思想將兩個字節的乘法運算進行拆分。

字節 y 與字節 a 的乘法運算過程可表示為，其中 y(a) 表示從 y 的乘法表中查詢與 x 相乘結果的操作：

y(a) = y * a

我們將字節 a 拆分成高4位（al) 與低 4 位 (ar) 兩個部分，即（其中 ⊕

為異或運算）：

a = (al << 4) ⊕ ar

這樣字節 a 就表示為 0-15 與 (0-15 << 4) 異或運算的結果了。于是原先的 y 與 a 的乘法運算可表示為：

y(a) = y(al << 4) ⊕ y(ar)

由于 ar 與 al 的范圍均為 0-15（0-1111），字節 y 與它們相乘的結果也就只有16個可能的值了 。這樣原先256 字節的字節 y 的乘法表就可以被 2 張 16 字節的乘法表替換了。

下面以根據本原多項式 x^8 + x^4 + x^3 + x^2 + 1 生成的 GF(2^8) 為例，分別通過查詢普通乘法表與使用拆分乘法表來演示 16 * 100 的計算過程。

16 的完整乘法表為：

table = [0 16 32 48 64 80 96 112 128 144 160 176 192 208 224 240
 29 13 61 45 93 77 125 109 157 141 189 173 221 205 253 237 58 42 26 10 122 
106 90 74 186 170 154 138 250 234 218 202 39 55 7 23 103 119 71 87 167 183
 135 151 231 247 199 215 116 100 84 68 52 36 20 4 244 228 212 196 180 164 
148 132 105 121 73 89 41 57 9 25 233 249 201 217 169 185 137 153 78 94 110 
126 14 30 46 62 206 222 238 254 142 158 174 190 83 67 115 99 19 3 51 35 211
 195 243 227 147 131 179 163 232 248 200 216 168 184 136 152 104 120 72 88 
40 56 8 24 245 229 213 197 181 165 149 133 117 101 85 69 53 37 21 5 210 194 
242 226 146 130 178 162 82 66 114 98 18 2 50 34 207 223 239 255 143 159 175 
191 79 95 111 127 15 31 47 63 156 140 188 172 220 204 252 236 28 12 60 44 
92 76 124 108 129 145 161 177 193 209 225 241 1 17 33 49 65 81 97 113 166 
182 134 150 230 246 198 214 38 54 6 22 102 118 70 86 187 171 155 139 251 
235 219 203 59 43 27 11 123 107 91 75]

計算 16 * 100 可以直接查表得到：

table[100] = 14

16 的低4位乘法表，也就是16 與 0-15 的乘法結果：

lowtable = [0 16 32 48 64 80 96 112 128 144 160 176 192 208 224 240]

16 的高4位乘法表，為16 與 0-15 << 4 的乘法結果：

hightable = [0 29 58 39 116 105 78 83 232 245 210 207 156 129 166 187]

將 100 （01100100）拆分，則：

100 = 0110 << 4 ⊕ 0100

在低位表中查詢 0100（4），得：

lowtable[4] = 64

在高位表中查詢0110 （6），得：

hightable[6] = 78

將兩個查詢結果異或：

result = 64 ^ 78 = 1000000 ^ 1001110 = 1110 = 14

從上面的對比中，我們不難發現采用SIMD的新算法提高查表速度主要表現在兩個方面：

1. 減少了乘法表大小；

2. 提高查表并行度（從1個字節到16甚至32個字節）

采用 SIMD 指令在大大降低了乘法表的規模的同時多了一次查表操作以及異或運算。由于新的乘法表每一部分只有 16 字節，我們可以順利的將其放置于 XMM 寄存器中，從而利用 SIMD 指令集提供的指令來進行數據向量運算，將原先的逐字節查表改進為并行的對 16 字節進行查表，同時異或操作也是 16 字節并行的。除此之外，由于乘法表的總體規模的下降，在編碼過程中的緩存污染也被大大減輕了，關于緩存的問題我們會在接下來的小節中進行更細致的分析。

以上的計算過程以單個字節作為例子，下面我們一同來分析利用 SIMD 技術對多個字節進行運算的過程。基本步驟如下：
拆分保存原始數據的 XMM 寄存器中的數據向量，分別存儲于不同的 XMM 寄存器中

根據拆分后的數據向量對乘法表進行重排，即得到查表結果。我們可以將乘法表理解為按順序排放的數組，數組長度為 16，查表的過程可以理解為將拆分后的數據（數據范圍為 0-15 ）作為索引對乘法表數組進行重新排序。這樣我們就可以通過排序指令完成查表操作了將重排后的結果進行異或，得到最終的運算結果

以下是偽代碼：

*// 將原始數據的右移4bit*

d2 = raw_data >> 4

*// 將右移后的數據的每字節與15（即1111）做AND操作,得到數據高位*

high_data = d2 AND 1111

*// 原始數據的每字節與15（即1111）做AND操作,得到數據低位*

low_data = raw_data AND 1111

*// 以數據作為索引對乘法表進行了重排*

for i, b = range low_data { low_ret[i]=low_table[b]}

for i, b = range high_data {high_ret[i]=high_table[b]}

*// 異或兩部分結果得到最終數據*

ret = low_ret XOR high_ret

需要注意的是，要使用 SIMD 加速有限域運算，對 CPU 的最低要求是支持 SSSE3 擴展指令集。另外為了充分提高效率，我們應該事先對數據進行內存對齊操作，在 SSSE3 下我們需要將數據對齊到 16 Bytes，否則我們只能使用非對齊指令進行數據的讀取和寫入。在這一點上比較特殊的是 Go 語言， 一方面 Go 支持直接調用匯編函數這為使用 SIMD 指令集提供了語言上的支持；但另外一方面 Golang 又隱藏了內存申請的細節，這使得指定內存對齊操作不可控，雖然我們也可以通過 cgo 或者匯編來實現，但這增加額外的負擔。所幸，對于 CPU 來說一個 Cache line 的大小為64byte，這在一定程度上可以幫助我們減少非對齊讀寫帶來的懲罰。另外，根據Golang 的內存對齊算法，對于較大的數據塊，Golang 是會自動對齊到 32 byte 的，因此對齊或非對齊指令的執行效果是一致的。

二 寫緩存友好代碼

緩存優化通過兩方面進行，其一是減少緩存污染；其二是提高緩存命中率。在嘗試做到這兩點之前，我們先來分析緩存的基本工作原理。

CPU 緩存的默認工作模式是 Write-Back, 即每一次讀寫內存數據都需要先寫入緩存。上文提到的 Cache line 即為緩存工作的基本單位，其大小為固定的 64 byte ，也就說哪怕從內存中讀取 1字節的數據，CPU 也會將其余的63 字節帶入緩存。這樣設計的原因主要是為了提高緩存的時間局域性，因為所要執行的數據大小通常遠遠超過這個數字，提前將數據讀取至緩存有利于接下來的數據在緩存中被命中。

2.1 矩陣運算分塊

矩陣運算的循環迭代中都用到了行與列，因此原始數據矩陣與編碼矩陣的訪問總有一方是非連續的，通過簡單的循環交換并不能改善運算的空間局域性。因此我們通過分塊的方法來提高時間局域性來減少緩存缺失。

分塊算法不是對一個數組的整行或整列進行操作，而是對其子矩陣進行操作，目的是在緩存中的數據被替換之前，最大限度的利用它。

分塊的尺寸不宜過大，太大的分塊無法被裝進緩存；另外也不能過小，太小的分塊導致外部邏輯的調用次數大大上升，產生了不必要的函數調用開銷，而且也不能充分利用緩存空間。

2.2 減少緩存污染

不難發現的是，編碼矩陣中的系數并不會完全覆蓋整個 GF(2^8)，例如 10+4 的編碼方案中，編碼矩陣中校驗矩陣大小為 4×10，編碼系數至多（可能會有重復）有10×4=40 個。因此我們可以事先進行一個乘法表初始化的過程，比如生成一個新的二維數組來存儲編碼系數的乘法表。縮小表的范圍可以在讀取表的過程中對緩存的污染。

另外在定義方法集時需要注意的是避免結構體中的元素浪費。避免將不必要的參數扔進結構體中，如果每一個方法僅使用其中若干個元素，則其他元素白白侵占了緩存空間。

三、 指令級并行與數據級并行的深入優化

本節主要介紹如何利用 AVX/AVX2 指令集以及指令級并行優化來進一步提高性能表現。除此之外，我們還可以對匯編代碼進行微調以取得微小的提升。比如，盡量避免使用 R8-R15 這 8 個寄存器，因為指令解碼會比其他通用寄存器多一個字節。但很多匯編優化細節是和 CPU 架構設計相關的，書本上甚至 Intel 提供的手冊也并不能提供最準確的指導（因為有滯后性），而且這些操作帶來的效益并不顯著，在這里就不做重點說明了。

3.1 利用 AVX2

在上文中我們已經知道如何將乘法表拆分成 128bits 的大小以適應 XMM 寄存器，那么對于 AVX 指令集來說，要充分發揮其作用，需要將乘法表復制到 256 bit 的 YMM 寄存器。為了做到這一點，我們可以利用 XMM 寄存器為 YMM 寄存器的低位這一特性，僅使用一條指令來完成表的復制（Intel 風格）：
vinserti128 ymm0, ymm0, xmm0, 1

這條指令作用是將 xmm0 寄存器中的數據拷貝到 ymm0 中,而剩余 128 位數據通過 ymm0 得到，其中立即數 1 表明 xmm0 拷貝的目的地是 ymm0 的高位。這條指令提供了兩個 source operand（源操作數）以及一個 destination operand（目標操作數），我們在這里使用 ymm0 寄存器同時作為源操作數
和目標操作數
來實現了表的復制操作。接下來我們便可以使用與 SSSE3 下同樣的方式來進行單指令 32 byte 的編碼運算過程了。

由于使用了 SSE 與 AVX 這兩種擴展指令集，我們需要避免 AVX-SSE Transition Penalties[3]。之所以會有這種性能懲罰主要是由于 SSE 指令對 YMM 寄存器的高位一無所知，SSE 指令與 AVX 指令的混用會導致機器不斷的執行 YMM 寄存器的高位保存與恢復，這大大影響了性能表現。如果對指令不熟悉，難以避免指令混用，那么可以在 RET 前使用 VZEROUPPER 指令來清空 YMM 寄存器的高位。

3.2 指令級并行 (ILP) 優化

程序分支指令的開銷并不僅僅為指令執行所需要的周期，因為它們可能影響前端流水線和內部緩存的內容。我們可以通過如下技巧來減少分支指令對性能的影響，并且提高分支預測單元的準確性：

1. 少的使用分支指令

2. 當貫穿 (fall-through) 更可能被執行時，使用向前條件跳轉

3. 當貫穿代碼不太可能被執行時，使用向后條件跳轉

向前跳轉經常用在檢查函數參數的代碼塊中，如果我們避免了傳入長度為 0 的數據切片，這樣可以在匯編中去掉相關的分支判斷。在我的代碼中僅有一條向后條件跳轉指令，用在循環代碼塊的底部。需要注意的是，以上 2 、 3 點中的優化方法是為了符合靜態分支預測算法的要求，然而在市場上基于硬件動態預測方法等處理器占主導地位，因此這兩點優化可能并不會起到提高分支預測準確度的作用，更多的是良好的編程習慣的問題。

對于 CPU 的執行引擎來說，其往往包含多個執行單元實例，這是執行引擎并發執行多個微操做的基本原理。另外 CPU 內核的調度器下會掛有多個端口，這意味著每個周期調度器可以給執行引擎分發多個微操作。因此我們可以利用循環展開來提高指令級并行的可能性。

循環展開就是將循環體復制多次，同時調整循環的終止代碼。由于它減少了分支判斷的次數，因此可以將來自不同迭代的指令放在一起調度。

當然，如果循環展開知識簡單地進行指令復制，最后使用的都是同一組寄存器，可能會妨礙對循環的有效調度。因此我們應當合理分配寄存器的使用。另外，如果循環規模較大，會導致指令緩存的缺失率上升。Intel 的優化手冊中指出，循環體不應當超過 500 條指令。[4]

四、 小結

以上內容較為完整的還原了糾刪碼引擎的實現過程，涉及到了較多的數學和硬件層面的知識，對于大部分工程師來說可能相對陌生，我們希望通過本系列文章的介紹能夠為大家的工程實踐提供些許幫助。但受限于篇幅，很多內容無法全面展開。比如，部分數學工具的理論與證明并沒有得到詳細的解釋，還需要讀者通過其他專業資料的來進行更深入的學習。

附錄：
Galois Fields and Cyclic Codes
http://user.xmission.com/~rimrock/Documents/Galois%20Fields%20and%20Cyclic%20Codes.pdf

有限域相關計算 https://github.com/templexxx/reedsolomon/tree/master/mathtools

Avoiding AVX-SSE Transition Penalties
https://software.intel.com/en-us/articles/avoiding-avx-sse-transition-penalties

Intel 64 and IA-32 Architectures Optimization Reference Manual :3.4.2.6 Optimization for Decoded ICache