Erasure coding implementation ============================= This document provides information about how [erasure code][1] has been implemented into ec translator. It describes the algorithm used and the optimizations made, but it doesn't contain a full description of the mathematical background needed to understand erasure coding in general. It either describes the other parts of ec not directly related to the encoding/decoding procedure, like synchronization or fop management. Introduction ------------ EC is based on [Reed-Solomon][2] erasure code. It's a very old code. It's not considered the best one nowadays, but is good enough and it's one of the few codes that is not covered by any patent and can be freely used. To define the Reed-Solomon code we use 3 parameters: * __Key fragments (K)__ It represents the minimum number of healthy fragments that will be needed to be able to recover the original data. Any subset of K out of the total number of fragments will serve. * __Redundancy fragments (R)__ It represents the number of extra fragments to compute for each original data block. This value determines how many fragments can be lost before being unable to recover the original data. * __Fragment size (S)__ This determines the size of each fragment. The original data block size is computed as S * K. Currently this values is fixed to 512 bytes. * __Total number of fragments (N = K + R)__ This isn't a real parameter but it will be useful to simplify the following descriptions. From the point of view of the implementation, it only consists on matrix multiplications. There are two kinds of matrices to use for Reed-Solomon: * __[Systematic][3]__ This kind of matrix has the particularity that K of the encoded fragments are simply a copy of the original data, divided into K pieces. Thus no real encoding needs to be done for them and only the R redundancy fragments need to be computed. This kind of matrices contain one KxK submatrix that is the [identity matrix][4]. * __Non-systematic__ This kind of matrix doesn't contain an identity submatrix. This means that all of the N fragments need to be encoded, requiring more computation. On the other hand, these matrices have some nice properties that allow faster implementations of some algorithms, like the matrix inversion used to decode the data. Another advantage of non-systematic matrices is that the decoding time is constant, independently of how many fragments are lost, while systematic approach can suffer from performance degradation when one fragment is lost. All non-systematic matrices can be converted to systematic ones, but then we lose the good properties of the non-systematic. We have to choose betwee best peek performance (systematic) and performance stability (non-systematic). Encoding procedure ------------------ To encode a block of data we need a KxN matrix where each subset of K rows is [linearly independent][5]. In other words, the determinant of each KxK submatrix is not 0. There are some known ways to obtain this kind of matrices. EC uses a small variation of a matrix known as [Vandermonde Matrix][6] where each element of the matrix is defined as: a(i, j) = i ^ (K - j) where i is the row from 1 to N, and j is the column from 1 to K. This is exactly the Vandermonde Matrix but with the elements of each row in reverse order. This change is made to be able to implement a small optimization in the matrix multiplication. Once we have the matrix, we only need to compute the multiplication of this matrix by a vector composed of K elements of data coming from the original data block. / \ / \ | 1 1 1 1 1 | / \ | a + b + c + d + e = t | | 16 8 4 2 1 | | a | | 16a + 8b + 4c + 2d + e = u | | 81 27 9 3 1 | | b | = | 81a + 27b + 9c + 3d + e = v | | 256 64 16 4 1 | * | c | | 256a + 64b + 16c + 4d + e = w | | 625 125 25 5 1 | | d | | 625a + 125b + 25c + 5d + e = x | | 1296 216 36 6 1 | | e | | 1296a + 216b + 36c + 6d + e = y | | 2401 343 49 7 1 | \ / | 2401a + 343b + 49c + 7d + e = z | \ / \ / The optimization that can be done here is this: 16a + 8b + 4c + 2d + e = 2(2(2(2a + b) + c) + d) + e So all the multiplications are always by the number of the row (2 in this case) and we don't need temporal storage for intermediate results: a *= 2 a += b a *= 2 a += c a *= 2 a += d a *= 2 a += e Once we have the result vector, each element is a fragment that needs to be stored in a separate place. Decoding procedure ------------------ To recover the data we need exactly K of the fragments. We need to know which K fragments we have (i.e. the original row number from which each fragment was calculated). Once we have this data we build a square KxK matrix composed by the rows corresponding to the given fragments and invert it. With the inverted matrix, we can recover the original data by multiplying it with the vector composed by the K fragments. In our previous example, if we consider that we have recovered fragments t, u, v, x and z, corresponding to rows 1, 2, 3, 5 and 7, we can build the following matrix: / \ | 1 1 1 1 1 | | 16 8 4 2 1 | | 81 27 9 3 1 | | 625 125 25 5 1 | | 2401 343 49 7 1 | \ / And invert it: / \ | 1/48 -1/15 1/16 -1/48 1/240 | | -17/48 16/15 -15/16 13/48 -11/240 | | 101/48 -86/15 73/16 -53/48 41/240 | | -247/48 176/15 -129/16 83/48 -61/240 | | 35/8 -7 35/8 -7/8 1/8 | \ / Multiplying it by the vector (t, u, v, x, z) we recover the original data (a, b, c, d, e): / \ / \ / \ | 1/48 -1/15 1/16 -1/48 1/240 | | t | | a | | -17/48 16/15 -15/16 13/48 -11/240 | | u | | b | | 101/48 -86/15 73/16 -53/48 41/240 | * | v | = | c | | -247/48 176/15 -129/16 83/48 -61/240 | | x | | d | | 35/8 -7 35/8 -7/8 1/8 | | z | | e | \ / \ / \ / Galois Field ------------ This encoding/decoding procedure is quite complex to compute using regular mathematical operations and it's not well suited for what we want to do (note that matrix elements can grow unboundly). To solve this problem, exactly the same procedure is done inside a [Galois Field][7] of characteristic 2, which is a finite field with some interesting properties that make it specially useful for fast operations using computers. There are two main differences when we use this specific Galois Field: * __All regular additions are replaced by bitwise xor's__ For todays computers it's not really faster to execute an xor compared to an addition, however replacing additions by xor's inside a multiplication has many advantages (we will make use of this to optimize the multiplication). Another consequence of this change is that additions and substractions are really the same xor operation. * __The elements of the matrix are bounded__ The field uses a modulus that keep all possible elements inside a delimited region, avoiding really big numbers and fixing the number of bits needed to represent each value. In the current implementation EC uses 8 bits per field element. It's very important to understand how multiplications are computed inside a Galois Field to be able to understand how has it been optimized. We'll start with a simple 'old school' multiplication but in base 2. For example, if we want to multiply 7 * 5 (111b * 101b in binary), we do the following: 1 1 1 (= 7) * 1 0 1 (= 5) ----------- 1 1 1 (= 7) + 0 0 0 (= 0) + 1 1 1 (= 7) ----------- 1 0 0 0 1 1 (= 35) This is quite simple. Note that the addition of the third column generates a carry that is propagated to all the other left columns. The next step is to define the modulus of the field. Suppose we use 11 as the modulus. Then we convert the result into an element of the field by dividing by the modulus and taking the residue. We also use the 'old school' method in binary: 1 0 0 0 1 1 (= 35) | 1 0 1 1 (= 11) - 0 0 0 0 ---------------- --------- 0 1 1 (= 3) 1 0 0 0 1 - 1 0 1 1 ----------- 0 0 1 1 0 1 - 1 0 1 1 ------------- 0 0 1 0 (= 2) So, 7 * 5 in a field with modulus 11 is 2. Note that the main objective in each iteration of the division is to make higher bits equal to 0 when possible (if it's not possible in one iteration, it will be zeroed on the next). If we do the same but changing additions with xors we get this: 1 1 1 (= 7) * 1 0 1 (= 5) ----------- 1 1 1 (= 7) x 0 0 0 (= 0) x 1 1 1 (= 7) ----------- 1 1 0 1 1 (= 27) In this case, the xor of the third column doesn't generate any carry. Now we need to divide by the modulus. We can also use 11 as the modulus since it still satisfies the needed conditions to work on a Galois Field of characteristic 2 with 3 bits: 1 1 0 1 1 (= 27) | 1 0 1 1 (= 11) x 1 0 1 1 ---------------- --------- 1 1 1 (= 7) 0 1 1 0 1 x 1 0 1 1 ----------- 0 1 1 0 1 x 1 0 1 1 ------------- 0 1 1 0 (= 6) Note that, in this case, to make zero the higher bit we need to consider the result of the xor operation, not the addition operation. So, 7 * 5 in a Galois Field of 3 bits with modulus 11 is 6. Optimization ------------ To compute all these operations in a fast way some methods have been traditionally used. Maybe the most common is the [lookup table][8]. The problem with this method is that it requires 3 lookups for each byte multiplication, greatly amplifying the needed memory bandwidth and making it difficult to take advantage of any SIMD support on the processor. What EC does to improve the performance is based on the following property (using the 3 bits Galois Field of the last example): A * B mod N = (A * b{2} * 4 mod N) x (A * b{1} * 2 mod N) x (A * b{0} mod N) This is basically a rewrite of the steps made in the previous example to multiply two numbers but moving the modulus calculation inside each intermediate result. What we can see here is that each term of the xor can be zeroed if the corresponding bit of B is 0, so we can ignore that factor. If the bit is 1, we need to compute A multiplied by a power of two and take the residue of the division by the modulus. We can precompute these values: A0 = A (we don't need to compute the modulus here) A1 = A0 * 2 mod N A2 = A1 * 2 mod N Having these values we only need to add those corresponding to bits set to 1 in B. Using our previous example: A = 1 1 1 (= 7) B = 1 0 1 (= 5) A0 = 1 1 1 (= 7) A1 = 1 1 1 * 1 0 mod 1 0 1 1 = 1 0 1 (= 5) A2 = 1 0 1 * 1 0 mod 1 0 1 1 = 0 0 1 (= 1) Since only bits 0 and 2 are 1 in B, we add A0 and A2: A0 + A2 = 1 1 1 x 0 0 1 = 1 1 0 (= 6) If we carefully look at what we are doing when computing each Ax, we see that we do two basic things: - Shift the original value one bit to the left - If the highest bit is 1, xor with the modulus Let's write this in a detailed way (representing each bit): Original value: a{2} a{1} a{0} Shift 1 bit: a{2} a{1} a{0} 0 If a{2} is 0 we already have the result: a{1} a{0} 0 If a{2} is 1 we need to xor with the modulus: 1 a{1} a{0} 0 x 1 0 1 1 = a{1} (a{0} x 1) 1 An important thing to see here is that if a{2} is 0, we can get the same result by xoring with all 0 instead of the modulus. For this reason we can rewrite the modulus as this: Modulus: a{2} 0 a{2} a{2} This means that the modulus will be 0 0 0 0 is a{2} is 0, so the value won't change, and it will be 1 0 1 1 if a{2} is 1, giving the correct result. So, the computation is simply: Original value: a{2} a{1} a{0} Shift 1 bit: a{2} a{1} a{0} 0 Apply modulus: a{1} (a{0} x a{2}) a{2} We can compute all Ax using this method. We'll get this: A0 = a{2} a{1} a{0} A1 = a{1} (a{0} x a{2}) a{2} A2 = (a{0} x a{2}) (a{1} x a{2}) a{1} Once we have all terms, we xor the ones corresponding to the bits set to 1 in B. In out example this will be A0 and A2: Result: (a{2} x a{0} x a{2}) (a{1} x a{1} x a{2}) (a{0} x a{1}) We can easily see that we can remove some redundant factors: Result: a{0} a{2} (a{0} x a{1}) This way we have come up with a simply set of equations to compute the multiplication of any number by 5. If A is 1 1 1 (= 7), the result must be 1 1 0 (= 6) using the equations, as we expected. If we try another numbe for A, like 0 1 0 (= 2), the result must be 0 0 1 (= 1). This seems a really fast way to compute the multiplication without using any table lookup. The problem is that this is only valid for B = 5. For other values of B another set of equations will be found. To solve this problem we can pregenerate the equations for all possible values of B. Since the Galois Field we use is small, this is feasible. One thing to be aware of is that, in general, two equations for different bits of the same B can share common subexpressions. This gives space for further optimizations to reduce the total number of xors used in the final equations for a given B. However this is not easy to find, since finding the smallest number of xors that give the correct result is an NP-Problem. For EC an exhaustive search has been made to find the best combinations for each possible value. Implementation -------------- All this seems great from the hardware point of view, but implementing this using normal processor instructions is not so easy because we would need a lot of shifts, ands, xors and ors to move the bits of each number to the correct position to compute the equation and then another shift to put each bit back to its final place. For example, to implement the functions to multiply by 5, we would need something like this: Bit 2: T2 = (A & 1) << 2 Bit 1: T1 = (A & 4) >> 1 Bit 0: T0 = ((A >> 1) x A) & 1 Result: T2 + T1 + T0 This doesn't look good. So here we make a change in the way we get and process the data: instead of reading full numbers into variables and operate with them afterwards, we use a single independent variable for each bit of the number. Assume that we can read and write independent bits from memory (later we'll see how to solve this problem when this is not possible). In this case, the code would look something like this: Bit 2: T2 = Mem[2] Bit 1: T1 = Mem[1] Bit 0: T0 = Mem[0] Computation: T1 ^= T0 Store result: Mem[2] = T0 Mem[1] = T2 Mem[0] = T1 Note that in this case we handle the final reordering of bits simply by storing the right variable to the right place, without any shifts, ands nor ors. In fact we only have memory loads, memory stores and xors. Note also that we can do all the computations directly using the variables themselves, without additional storage. This true for most of the values, but in some cases an additional one or two temporal variables will be needed to store intermediate results. The drawback of this approach is that additions, that are simply a xor of two numbers will need as many xors as bits are in each number. SIMD optimization ----------------- So we have a good way to compute the multiplications, but even using this we'll need several operations for each byte of the original data. We can improve this by doing multiple multiplications using the same set of instructions. With the approach taken in the implementation section, we can see that in fact it's really easy to add SIMD support to this method. We only need to store in each variable one bit from multiple numbers. For example, when we load T2 from memory, instead of reading the bit 2 of the first number, we can read the bit 2 of the first, second, third, fourth, ... numbers. The same can be done when loading T1 and T0. Obviously this needs to have a special encoding of the numbers into memory to be able to do that in a single operation, but since we can choose whatever encoding we want for EC, we have chosen to have exactly that. We interpret the original data as a stream of bits, and we split it into subsequences of length L, each containing one bit of a number. Every S subsequences form a set of numbers of S bits that are encoded and decoded as a single group. This repeats for any remaining data. For example, in a simple case with L = 8 and S = 3, the original data would contain something like this (interpreted as a sequence of bits, offsets are also bit-based): Offset 0: a{0} b{0} c{0} d{0} e{0} f{0} g{0} h{0} Offset 8: a{1} b{1} c{1} d{1} e{1} f{1} g{1} h{1} Offset 16: a{2} b{2} c{2} d{2} e{2} f{2} g{2} h{2} Offset 24: i{0} j{0} k{0} l{0} m{0} n{0} o{0} p{0} Offset 32: i{1} j{1} k{1} l{1} m{1} n{1} o{1} p{1} Offset 40: i{2} j{2} k{2} l{2} m{2} n{2} o{2} p{2} Note: If the input file is not a multiple of S * L, 0-padding is done. Here we have 16 numbers encoded, from A to P. This way we can easily see that reading the first byte of the file will read all bits 0 of number A, B, C, D, E, F, G and H. The same happens with bits 1 and 2 when we read the second and third bytes respectively. Using this encoding and the implementation described above, we can see that the same set of instructions will be computing the multiplication of 8 numbers at the same time. This can be further improved if we use L = 64 with 64 bits variables on 64-bits processor. It's even faster if we use L = 128 using SSE registers or L = 256 using AVX registers on Intel processors. Currently EC uses L = 512 and S = 8. This means that numbers are packed in blocks of 512 bytes and gives space for even bigger processor registers up to 512 bits. Conclusions ----------- This method requires a single variable/processor register for each bit. This can be challenging if we want to avoid additional memory accesses, even if we use modern processors that have many registers. However, the implementation we chose for the Vandermonde Matrix doesn't require temporary storage, so we don't need a full set of 8 new registers (one for each bit) to store partial computations. Additionally, the computation of the multiplications requires, at most, 2 extra registers, but this is afordable. Xors are a really fast operation in modern processors. Intel CPU's can dispatch up to 3 xors per CPU cycle if there are no dependencies with ongoing previous instructions. Worst case is 1 xor per cycle. So, in some configurations, this method could be very near to the memory speed. Another interesting thing of this method is that all data it needs to operate is packed in small sequential blocks of memory, meaning that it can take advantage of the faster internal CPU caches. Results ------- For the particular case of 8 bits, EC can compute each multiplication using 12.8 xors on average (without counting 0 and 1 that do not require any xor). Some numbers require less, like 2 that only requires 3 xors. Having all this, we can check some numbers to see the performance of this method. Maybe the most interesting thing is the average number of xors needed to encode a single byte of data. To compute this we'll need to define some variables: * K: Number of data fragments * R: Number of redundancy fragments * N: K + R * B: Number of bits per number * A: Average number of xors per number * Z: Bits per CPU register (can be up to 256 for AVX registers) * X: Average number of xors per CPU cycle * L: Average cycles per load * S: Average cycles per store * G: Core speed in Hz _Total number of bytes processed for a single matrix multiplication_: * __Read__: K * B * Z / 8 * __Written__: N * B * Z / 8 _Total number of memory accesses_: * __Loads__: K * B * N * __Stores__: B * N > We need to read the same K * B * Z bits, in registers of Z bits, N > times, one for each row of the matrix. However the last N - 1 reads > could be made from the internal CPU caches if conditions are good. _Total number of operations_: * __Additions__: (K - 1) * N * __Multiplications__: K * N __Total number of xors__: B * (K - 1) * N + A * K * N = N * ((A + B) * K - B) __Xors per byte__: 8 * N * ((A + B) * K - B) / (K * B * Z) __CPU cycles per byte__: 8 * N * ((A + B) * K - B) / (K * B * Z * X) + 8 * L * N / Z + (loads) 8 * S * N / (K * Z) (stores) __Bytes per second__: G / {CPU cycles per byte} Some xors per byte numbers for specific configurations (B=8): Z=64 Z=128 Z=256 K=2/R=1 0.79 0.39 0.20 K=4/R=2 1.76 0.88 0.44 K=4/R=3 2.06 1.03 0.51 K=8/R=3 3.40 1.70 0.85 K=8/R=4 3.71 1.86 0.93 K=16/R=4 6.34 3.17 1.59 [1]: https://en.wikipedia.org/wiki/Erasure_code [2]: https://en.wikipedia.org/wiki/Reed%E2%80%93Solomon_error_correction [3]: https://en.wikipedia.org/wiki/Systematic_code [4]: https://en.wikipedia.org/wiki/Identity_matrix [5]: https://en.wikipedia.org/wiki/Linear_independence [6]: https://en.wikipedia.org/wiki/Vandermonde_matrix [7]: https://en.wikipedia.org/wiki/Finite_field [8]: https://en.wikipedia.org/wiki/Finite_field_arithmetic#Implementation_tricks