1 /**********************************************************************
2 Copyright(c) 2011-2018 Intel Corporation All rights reserved.
3
4 Redistribution and use in source and binary forms, with or without
5 modification, are permitted provided that the following conditions
6 are met:
7 * Redistributions of source code must retain the above copyright
8 notice, this list of conditions and the following disclaimer.
9 * Redistributions in binary form must reproduce the above copyright
10 notice, this list of conditions and the following disclaimer in
11 the documentation and/or other materials provided with the
12 distribution.
13 * Neither the name of Intel Corporation nor the names of its
14 contributors may be used to endorse or promote products derived
15 from this software without specific prior written permission.
16
17 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
18 "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
19 LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
20 A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
21 OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
22 SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
23 LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
24 DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
25 THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
26 (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
27 OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
28 **********************************************************************/
29
30 #include <stdio.h>
31 #include <stdlib.h>
32 #include <string.h>
33 #include <getopt.h>
34 #include "erasure_code.h" // use <isa-l.h> instead when linking against installed
35 #include "test.h"
36
37 #define MMAX 255
38 #define KMAX 255
39
40 typedef unsigned char u8;
41 int verbose = 0;
42
43 int
usage(void)44 usage(void)
45 {
46 fprintf(stderr,
47 "Usage: ec_piggyback_example [options]\n"
48 " -h Help\n"
49 " -k <val> Number of source fragments\n"
50 " -p <val> Number of parity fragments\n"
51 " -l <val> Length of fragments\n"
52 " -e <val> Simulate erasure on frag index val. Zero based. Can be repeated.\n"
53 " -v Verbose\n"
54 " -b Run timed benchmark\n"
55 " -s Toggle use of sparse matrix opt\n"
56 " -r <seed> Pick random (k, p) with seed\n");
57 exit(0);
58 }
59
60 // Cauchy-based matrix
61 void
gf_gen_full_pb_cauchy_matrix(u8 * a,int m,int k)62 gf_gen_full_pb_cauchy_matrix(u8 *a, int m, int k)
63 {
64 int i, j, p = m - k;
65
66 // Identity matrix in top k x k to indicate a symmetric code
67 memset(a, 0, k * m);
68 for (i = 0; i < k; i++)
69 a[k * i + i] = 1;
70
71 for (i = k; i < (k + p / 2); i++) {
72 for (j = 0; j < k / 2; j++)
73 a[k * i + j] = gf_inv(i ^ j);
74 for (; j < k; j++)
75 a[k * i + j] = 0;
76 }
77 for (; i < m; i++) {
78 for (j = 0; j < k / 2; j++)
79 a[k * i + j] = 0;
80 for (; j < k; j++)
81 a[k * i + j] = gf_inv((i - p / 2) ^ (j - k / 2));
82 }
83
84 // Fill in mixture of B parity depending on a few localized A sources
85 int r = 0, c = 0;
86 int repeat_len = k / (p - 2);
87 int parity_rows = p / 2;
88
89 for (i = 1 + k + parity_rows; i < m; i++, r++) {
90 if (r == (parity_rows - 1) - ((k / 2 % (parity_rows - 1))))
91 repeat_len++;
92
93 for (j = 0; j < repeat_len; j++, c++)
94 a[k * i + c] = gf_inv((k + 1) ^ c);
95 }
96 }
97
98 // Vandermonde based matrix - not recommended due to limits when invertable
99 void
gf_gen_full_pb_vand_matrix(u8 * a,int m,int k)100 gf_gen_full_pb_vand_matrix(u8 *a, int m, int k)
101 {
102 int i, j, p = m - k;
103 unsigned char q, gen = 1;
104
105 // Identity matrix in top k x k to indicate a symmetric code
106 memset(a, 0, k * m);
107 for (i = 0; i < k; i++)
108 a[k * i + i] = 1;
109
110 for (i = k; i < (k + (p / 2)); i++) {
111 q = 1;
112 for (j = 0; j < k / 2; j++) {
113 a[k * i + j] = q;
114 q = gf_mul(q, gen);
115 }
116 for (; j < k; j++)
117 a[k * i + j] = 0;
118 gen = gf_mul(gen, 2);
119 }
120 gen = 1;
121 for (; i < m; i++) {
122 q = 1;
123 for (j = 0; j < k / 2; j++) {
124 a[k * i + j] = 0;
125 }
126 for (; j < k; j++) {
127 a[k * i + j] = q;
128 q = gf_mul(q, gen);
129 }
130 gen = gf_mul(gen, 2);
131 }
132
133 // Fill in mixture of B parity depending on a few localized A sources
134 int r = 0, c = 0;
135 int repeat_len = k / (p - 2);
136 int parity_rows = p / 2;
137
138 for (i = 1 + k + parity_rows; i < m; i++, r++) {
139 if (r == (parity_rows - 1) - ((k / 2 % (parity_rows - 1))))
140 repeat_len++;
141
142 for (j = 0; j < repeat_len; j++)
143 a[k * i + c++] = 1;
144 }
145 }
146
147 void
print_matrix(int m,int k,unsigned char * s,const char * msg)148 print_matrix(int m, int k, unsigned char *s, const char *msg)
149 {
150 int i, j;
151
152 printf("%s:\n", msg);
153 for (i = 0; i < m; i++) {
154 printf("%3d- ", i);
155 for (j = 0; j < k; j++) {
156 printf(" %2x", 0xff & s[j + (i * k)]);
157 }
158 printf("\n");
159 }
160 printf("\n");
161 }
162
163 void
print_list(int n,unsigned char * s,const char * msg)164 print_list(int n, unsigned char *s, const char *msg)
165 {
166 int i;
167 if (!verbose)
168 return;
169
170 printf("%s: ", msg);
171 for (i = 0; i < n; i++)
172 printf(" %d", s[i]);
173 printf("\n");
174 }
175
176 static int
177 gf_gen_decode_matrix(u8 *encode_matrix, u8 *decode_matrix, u8 *invert_matrix, u8 *temp_matrix,
178 u8 *decode_index, u8 *frag_err_list, int nerrs, int k, int m);
179
180 int
main(int argc,char * argv[])181 main(int argc, char *argv[])
182 {
183 int i, j, m, c, e, ret;
184 int k = 10, p = 4, len = 8 * 1024; // Default params
185 int nerrs = 0;
186 int benchmark = 0;
187 int sparse_matrix_opt = 1;
188
189 // Fragment buffer pointers
190 u8 *frag_ptrs[MMAX];
191 u8 *parity_ptrs[KMAX];
192 u8 *recover_srcs[KMAX];
193 u8 *recover_outp[KMAX];
194 u8 frag_err_list[MMAX];
195
196 // Coefficient matrices
197 u8 *encode_matrix, *decode_matrix;
198 u8 *invert_matrix, *temp_matrix;
199 u8 *g_tbls;
200 u8 decode_index[MMAX];
201
202 if (argc == 1)
203 for (i = 0; i < p; i++)
204 frag_err_list[nerrs++] = rand() % (k + p);
205
206 while ((c = getopt(argc, argv, "k:p:l:e:r:hvbs")) != -1) {
207 switch (c) {
208 case 'k':
209 k = atoi(optarg);
210 break;
211 case 'p':
212 p = atoi(optarg);
213 break;
214 case 'l':
215 len = atoi(optarg);
216 if (len < 0)
217 usage();
218 break;
219 case 'e':
220 e = atoi(optarg);
221 frag_err_list[nerrs++] = e;
222 break;
223 case 'r':
224 srand(atoi(optarg));
225 k = (rand() % MMAX) / 4;
226 k = (k < 2) ? 2 : k;
227 p = (rand() % (MMAX - k)) / 4;
228 p = (p < 2) ? 2 : p;
229 for (i = 0; i < k && nerrs < p; i++)
230 if (rand() & 1)
231 frag_err_list[nerrs++] = i;
232 break;
233 case 'v':
234 verbose++;
235 break;
236 case 'b':
237 benchmark = 1;
238 break;
239 case 's':
240 sparse_matrix_opt = !sparse_matrix_opt;
241 break;
242 case 'h':
243 default:
244 usage();
245 break;
246 }
247 }
248 m = k + p;
249
250 // Check for valid parameters
251 if (m > (MMAX / 2) || k > (KMAX / 2) || m < 0 || p < 2 || k < 1) {
252 printf(" Input test parameter error m=%d, k=%d, p=%d, erasures=%d\n", m, k, p,
253 nerrs);
254 usage();
255 }
256 if (nerrs > p) {
257 printf(" Number of erasures chosen exceeds power of code erasures=%d p=%d\n", nerrs,
258 p);
259 }
260 for (i = 0; i < nerrs; i++) {
261 if (frag_err_list[i] >= m)
262 printf(" fragment %d not in range\n", frag_err_list[i]);
263 }
264
265 printf("ec_piggyback_example:\n");
266
267 /*
268 * One simple way to implement piggyback codes is to keep a 2x wide matrix
269 * that covers the how each parity is related to both A and B sources. This
270 * keeps it easy to generalize in parameters m,k and the resulting sparse
271 * matrix multiplication can be optimized by pre-removal of zero items.
272 */
273
274 int k2 = 2 * k;
275 int p2 = 2 * p;
276 int m2 = k2 + p2;
277 int nerrs2 = nerrs;
278
279 encode_matrix = malloc(m2 * k2);
280 decode_matrix = malloc(m2 * k2);
281 invert_matrix = malloc(m2 * k2);
282 temp_matrix = malloc(m2 * k2);
283 g_tbls = malloc(k2 * p2 * 32);
284
285 if (encode_matrix == NULL || decode_matrix == NULL || invert_matrix == NULL ||
286 temp_matrix == NULL || g_tbls == NULL) {
287 printf("Test failure! Error with malloc\n");
288 return -1;
289 }
290 // Allocate the src fragments
291 for (i = 0; i < k; i++) {
292 if (NULL == (frag_ptrs[i] = malloc(len))) {
293 printf("alloc error: Fail\n");
294 return -1;
295 }
296 }
297 // Allocate the parity fragments
298 for (i = 0; i < p2; i++) {
299 if (NULL == (parity_ptrs[i] = malloc(len / 2))) {
300 printf("alloc error: Fail\n");
301 return -1;
302 }
303 }
304
305 // Allocate buffers for recovered data
306 for (i = 0; i < p2; i++) {
307 if (NULL == (recover_outp[i] = malloc(len / 2))) {
308 printf("alloc error: Fail\n");
309 return -1;
310 }
311 }
312
313 // Fill sources with random data
314 for (i = 0; i < k; i++)
315 for (j = 0; j < len; j++)
316 frag_ptrs[i][j] = rand();
317
318 printf(" encode (m,k,p)=(%d,%d,%d) len=%d\n", m, k, p, len);
319
320 // Pick an encode matrix.
321 gf_gen_full_pb_cauchy_matrix(encode_matrix, m2, k2);
322
323 if (verbose)
324 print_matrix(m2, k2, encode_matrix, "encode matrix");
325
326 // Initialize g_tbls from encode matrix
327 ec_init_tables(k2, p2, &encode_matrix[k2 * k2], g_tbls);
328
329 // Fold A and B into single list of fragments
330 for (i = 0; i < k; i++)
331 frag_ptrs[i + k] = &frag_ptrs[i][len / 2];
332
333 if (!sparse_matrix_opt) {
334 // Standard encode using no assumptions on the encode matrix
335
336 // Generate EC parity blocks from sources
337 ec_encode_data(len / 2, k2, p2, g_tbls, frag_ptrs, parity_ptrs);
338
339 if (benchmark) {
340 struct perf start;
341 BENCHMARK(&start, BENCHMARK_TIME,
342 ec_encode_data(len / 2, k2, p2, g_tbls, frag_ptrs, parity_ptrs));
343 printf("ec_piggyback_encode_std: ");
344 perf_print(start, m2 * len / 2);
345 }
346 } else {
347 // Sparse matrix optimization - use fact that input matrix is sparse
348
349 // Keep an encode matrix with some zero elements removed
350 u8 *encode_matrix_faster, *g_tbls_faster;
351 encode_matrix_faster = malloc(m * k);
352 g_tbls_faster = malloc(k * p * 32);
353 if (encode_matrix_faster == NULL || g_tbls_faster == NULL) {
354 printf("Test failure! Error with malloc\n");
355 return -1;
356 }
357
358 /*
359 * Pack with only the part that we know are non-zero. Alternatively
360 * we could search and keep track of non-zero elements but for
361 * simplicity we just skip the lower quadrant.
362 */
363 for (i = k, j = k2; i < m; i++, j++)
364 memcpy(&encode_matrix_faster[k * i], &encode_matrix[k2 * j], k);
365
366 if (verbose) {
367 print_matrix(p, k, &encode_matrix_faster[k * k], "encode via sparse-opt");
368 print_matrix(p2 / 2, k2, &encode_matrix[(k2 + p2 / 2) * k2],
369 "encode via sparse-opt");
370 }
371 // Initialize g_tbls from encode matrix
372 ec_init_tables(k, p, &encode_matrix_faster[k * k], g_tbls_faster);
373
374 // Generate EC parity blocks from sources
375 ec_encode_data(len / 2, k, p, g_tbls_faster, frag_ptrs, parity_ptrs);
376 ec_encode_data(len / 2, k2, p, &g_tbls[k2 * p * 32], frag_ptrs, &parity_ptrs[p]);
377
378 if (benchmark) {
379 struct perf start;
380 BENCHMARK(&start, BENCHMARK_TIME,
381 ec_encode_data(len / 2, k, p, g_tbls_faster, frag_ptrs,
382 parity_ptrs);
383 ec_encode_data(len / 2, k2, p, &g_tbls[k2 * p * 32], frag_ptrs,
384 &parity_ptrs[p]));
385 printf("ec_piggyback_encode_sparse: ");
386 perf_print(start, m2 * len / 2);
387 }
388 }
389
390 if (nerrs <= 0)
391 return 0;
392
393 printf(" recover %d fragments\n", nerrs);
394
395 // Set frag pointers to correspond to parity
396 for (i = k2; i < m2; i++)
397 frag_ptrs[i] = parity_ptrs[i - k2];
398
399 print_list(nerrs2, frag_err_list, " frag err list");
400
401 // Find a decode matrix to regenerate all erasures from remaining frags
402 ret = gf_gen_decode_matrix(encode_matrix, decode_matrix, invert_matrix, temp_matrix,
403 decode_index, frag_err_list, nerrs2, k2, m2);
404
405 if (ret != 0) {
406 printf("Fail on generate decode matrix\n");
407 return -1;
408 }
409 // Pack recovery array pointers as list of valid fragments
410 for (i = 0; i < k2; i++)
411 if (decode_index[i] < k2)
412 recover_srcs[i] = frag_ptrs[decode_index[i]];
413 else
414 recover_srcs[i] = parity_ptrs[decode_index[i] - k2];
415
416 print_list(k2, decode_index, " decode index");
417
418 // Recover data
419 ec_init_tables(k2, nerrs2, decode_matrix, g_tbls);
420 ec_encode_data(len / 2, k2, nerrs2, g_tbls, recover_srcs, recover_outp);
421
422 if (benchmark) {
423 struct perf start;
424 BENCHMARK(&start, BENCHMARK_TIME,
425 ec_encode_data(len / 2, k2, nerrs2, g_tbls, recover_srcs, recover_outp));
426 printf("ec_piggyback_decode: ");
427 perf_print(start, (k2 + nerrs2) * len / 2);
428 }
429 // Check that recovered buffers are the same as original
430 printf(" check recovery of block {");
431 for (i = 0; i < nerrs2; i++) {
432 printf(" %d", frag_err_list[i]);
433 if (memcmp(recover_outp[i], frag_ptrs[frag_err_list[i]], len / 2)) {
434 printf(" Fail erasure recovery %d, frag %d\n", i, frag_err_list[i]);
435 return -1;
436 }
437 }
438 printf(" } done all: Pass\n");
439
440 return 0;
441 }
442
443 // Generate decode matrix from encode matrix and erasure list
444
445 static int
gf_gen_decode_matrix(u8 * encode_matrix,u8 * decode_matrix,u8 * invert_matrix,u8 * temp_matrix,u8 * decode_index,u8 * frag_err_list,int nerrs,int k,int m)446 gf_gen_decode_matrix(u8 *encode_matrix, u8 *decode_matrix, u8 *invert_matrix, u8 *temp_matrix,
447 u8 *decode_index, u8 *frag_err_list, int nerrs, int k, int m)
448 {
449 int i, j, p, r;
450 int nsrcerrs = 0;
451 u8 s, *b = temp_matrix;
452 u8 frag_in_err[MMAX];
453
454 memset(frag_in_err, 0, sizeof(frag_in_err));
455
456 // Order the fragments in erasure for easier sorting
457 for (i = 0; i < nerrs; i++) {
458 if (frag_err_list[i] < k)
459 nsrcerrs++;
460 frag_in_err[frag_err_list[i]] = 1;
461 }
462
463 // Construct b (matrix that encoded remaining frags) by removing erased rows
464 for (i = 0, r = 0; i < k; i++, r++) {
465 while (frag_in_err[r])
466 r++;
467 for (j = 0; j < k; j++)
468 b[k * i + j] = encode_matrix[k * r + j];
469 decode_index[i] = r;
470 }
471 if (verbose > 1)
472 print_matrix(k, k, b, "matrix to invert");
473
474 // Invert matrix to get recovery matrix
475 if (gf_invert_matrix(b, invert_matrix, k) < 0)
476 return -1;
477
478 if (verbose > 2)
479 print_matrix(k, k, invert_matrix, "matrix inverted");
480
481 // Get decode matrix with only wanted recovery rows
482 for (i = 0; i < nsrcerrs; i++) {
483 for (j = 0; j < k; j++) {
484 decode_matrix[k * i + j] = invert_matrix[k * frag_err_list[i] + j];
485 }
486 }
487
488 // For non-src (parity) erasures need to multiply encode matrix * invert
489 for (p = nsrcerrs; p < nerrs; p++) {
490 for (i = 0; i < k; i++) {
491 s = 0;
492 for (j = 0; j < k; j++)
493 s ^= gf_mul(invert_matrix[j * k + i],
494 encode_matrix[k * frag_err_list[p] + j]);
495
496 decode_matrix[k * p + i] = s;
497 }
498 }
499 if (verbose > 1)
500 print_matrix(nerrs, k, decode_matrix, "decode matrix");
501 return 0;
502 }
503