ceph/src/erasure-code/jerasure/jerasure/src/reed_sol.c

   1 /* *
   2  * Copyright (c) 2014, James S. Plank and Kevin Greenan
   3  * All rights reserved.
   4  *
   5  * Jerasure - A C/C++ Library for a Variety of Reed-Solomon and RAID-6 Erasure
   6  * Coding Techniques
   7  *
   8  * Revision 2.0: Galois Field backend now links to GF-Complete
   9  *
  10  * Redistribution and use in source and binary forms, with or without
  11  * modification, are permitted provided that the following conditions
  12  * are met:
  13  *
  14  *  - Redistributions of source code must retain the above copyright
  15  *    notice, this list of conditions and the following disclaimer.
  16  *
  17  *  - Redistributions in binary form must reproduce the above copyright
  18  *    notice, this list of conditions and the following disclaimer in
  19  *    the documentation and/or other materials provided with the
  20  *    distribution.
  21  *
  22  *  - Neither the name of the University of Tennessee nor the names of its
  23  *    contributors may be used to endorse or promote products derived
  24  *    from this software without specific prior written permission.
  25  *
  26  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
  27  * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
  28  * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
  29  * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
  30  * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
  31  * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
  32  * BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS
  33  * OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
  34  * AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
  35  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY
  36  * WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
  37  * POSSIBILITY OF SUCH DAMAGE.
  38  */
  39
  40 /* Jerasure's authors:
  41
  42    Revision 2.x - 2014: James S. Plank and Kevin M. Greenan
  43    Revision 1.2 - 2008: James S. Plank, Scott Simmerman and Catherine D. Schuman.
  44    Revision 1.0 - 2007: James S. Plank
  45  */
  46
  47 #include <stdio.h>
  48 #include <stdlib.h>
  49 #include <string.h>
  50 #include <assert.h>
  51
  52 #include <gf_complete.h>
  53 #include "galois.h"
  54 #include "jerasure.h"
  55 #include "reed_sol.h"
  56
  57 #define talloc(type, num) (type *) malloc(sizeof(type)*(num))
  58
  59 int *reed_sol_r6_coding_matrix(int k, int w)
  60 {
  61   int *matrix;
  62   int i, tmp;
  63
  64   if (w != 8 && w != 16 && w != 32) return NULL;
  65
  66   matrix = talloc(int, 2*k);
  67   if (matrix == NULL) return NULL;
  68
  69   for (i = 0; i < k; i++) matrix[i] = 1;
  70   matrix[k] = 1;
  71   tmp = 1;
  72   for (i = 1; i < k; i++) {
  73     tmp = galois_single_multiply(tmp, 2, w);
  74     matrix[k+i] = tmp;
  75   }
  76   return matrix;
  77 }
  78
  79 int *reed_sol_vandermonde_coding_matrix(int k, int m, int w)
  80 {
  81   int i, j;
  82   int *vdm, *dist;
  83
  84   vdm = reed_sol_big_vandermonde_distribution_matrix(k+m, k, w);
  85   if (vdm == NULL) return NULL;
  86   dist = talloc(int, m*k);
  87   if (dist == NULL) {
  88     free(vdm);
  89     return NULL;
  90   }
  91
  92   i = k*k;
  93   for (j = 0; j < m*k; j++) {
  94     dist[j] = vdm[i];
  95     i++;
  96   }
  97   free(vdm);
  98   return dist;
  99 }
 100
 101 static int prim08 = -1;
 102 static gf_t GF08;
 103
 104 void reed_sol_galois_w08_region_multby_2(char *region, int nbytes)
 105 {
 106   if (prim08 == -1) {
 107     prim08 = galois_single_multiply((1 << 7), 2, 8);
 108     if (!gf_init_hard(&GF08, 8, GF_MULT_BYTWO_b, GF_REGION_DEFAULT, GF_DIVIDE_DEFAULT,
 109                       prim08, 0, 0, NULL, NULL)) {
 110       fprintf(stderr, "Error: Can't initialize the GF for reed_sol_galois_w08_region_multby_2\n");
 111       assert(0);
 112     }
 113   }
 114   GF08.multiply_region.w32(&GF08, region, region, 2, nbytes, 0);
 115 }
 116
 117 static int prim16 = -1;
 118 static gf_t GF16;
 119
 120 void reed_sol_galois_w16_region_multby_2(char *region, int nbytes)
 121 {
 122   if (prim16 == -1) {
 123     prim16 = galois_single_multiply((1 << 15), 2, 16);
 124     if (!gf_init_hard(&GF16, 16, GF_MULT_BYTWO_b, GF_REGION_DEFAULT, GF_DIVIDE_DEFAULT,
 125                       prim16, 0, 0, NULL, NULL)) {
 126       fprintf(stderr, "Error: Can't initialize the GF for reed_sol_galois_w16_region_multby_2\n");
 127       assert(0);
 128     }
 129   }
 130   GF16.multiply_region.w32(&GF16, region, region, 2, nbytes, 0);
 131 }
 132
 133 static int prim32 = -1;
 134 static gf_t GF32;
 135
 136 void reed_sol_galois_w32_region_multby_2(char *region, int nbytes)
 137 {
 138   if (prim32 == -1) {
 139     prim32 = galois_single_multiply(((gf_val_32_t)1 << 31), 2, 32);
 140     if (!gf_init_hard(&GF32, 32, GF_MULT_BYTWO_b, GF_REGION_DEFAULT, GF_DIVIDE_DEFAULT,
 141                       prim32, 0, 0, NULL, NULL)) {
 142       fprintf(stderr, "Error: Can't initialize the GF for reed_sol_galois_w32_region_multby_2\n");
 143       assert(0);
 144     }
 145   }
 146   GF32.multiply_region.w32(&GF32, region, region, 2, nbytes, 0);
 147 }
 148
 149 int reed_sol_r6_encode(int k, int w, char **data_ptrs, char **coding_ptrs, int size)
 150 {
 151   int i;
 152
 153   /* First, put the XOR into coding region 0 */
 154
 155   memcpy(coding_ptrs[0], data_ptrs[0], size);
 156
 157   for (i = 1; i < k; i++) galois_region_xor(data_ptrs[i], coding_ptrs[0], size);
 158
 159   /* Next, put the sum of (2^j)*Dj into coding region 1 */
 160
 161   memcpy(coding_ptrs[1], data_ptrs[k-1], size);
 162
 163   for (i = k-2; i >= 0; i--) {
 164     switch (w) {
 165       case 8:  reed_sol_galois_w08_region_multby_2(coding_ptrs[1], size); break;
 166       case 16: reed_sol_galois_w16_region_multby_2(coding_ptrs[1], size); break;
 167       case 32: reed_sol_galois_w32_region_multby_2(coding_ptrs[1], size); break;
 168       default: return 0;
 169     }
 170
 171     galois_region_xor(data_ptrs[i], coding_ptrs[1], size);
 172   }
 173   return 1;
 174 }
 175
 176 int *reed_sol_extended_vandermonde_matrix(int rows, int cols, int w)
 177 {
 178   int *vdm;
 179   int i, j, k;
 180
 181   if (w < 30 && (1 << w) < rows) return NULL;
 182   if (w < 30 && (1 << w) < cols) return NULL;
 183
 184   vdm = talloc(int, rows*cols);
 185   if (vdm == NULL) { return NULL; }
 186
 187   vdm[0] = 1;
 188   for (j = 1; j < cols; j++) vdm[j] = 0;
 189   if (rows == 1) return vdm;
 190
 191   i=(rows-1)*cols;
 192   for (j = 0; j < cols-1; j++) vdm[i+j] = 0;
 193   vdm[i+j] = 1;
 194   if (rows == 2) return vdm;
 195
 196   for (i = 1; i < rows-1; i++) {
 197     k = 1;
 198     for (j = 0; j < cols; j++) {
 199       vdm[i*cols+j] = k;
 200       k = galois_single_multiply(k, i, w);
 201     }
 202   }
 203   return vdm;
 204 }
 205
 206 int *reed_sol_big_vandermonde_distribution_matrix(int rows, int cols, int w)
 207 {
 208   int *dist;
 209   int i, j, k;
 210   int sindex, srindex, siindex, tmp;
 211
 212   if (cols >= rows) return NULL;
 213
 214   dist = reed_sol_extended_vandermonde_matrix(rows, cols, w);
 215   if (dist == NULL) return NULL;
 216
 217   sindex = 0;
 218   for (i = 1; i < cols; i++) {
 219     sindex += cols;
 220
 221     /* Find an appropriate row -- where i,i != 0 */
 222     srindex = sindex+i;
 223     for (j = i; j < rows && dist[srindex] == 0; j++) srindex += cols;
 224     if (j >= rows) {   /* This should never happen if rows/w are correct */
 225       fprintf(stderr, "reed_sol_big_vandermonde_distribution_matrix(%d,%d,%d) - couldn't make matrix\n",
 226              rows, cols, w);
 227       assert(0);
 228     }
 229
 230     /* If necessary, swap rows */
 231     if (j != i) {
 232       srindex -= i;
 233       for (k = 0; k < cols; k++) {
 234         tmp = dist[srindex+k];
 235         dist[srindex+k] = dist[sindex+k];
 236         dist[sindex+k] = tmp;
 237       }
 238     }
 239
 240     /* If Element i,i is not equal to 1, multiply the column by 1/i */
 241
 242     if (dist[sindex+i] != 1) {
 243       tmp = galois_single_divide(1, dist[sindex+i], w);
 244       srindex = i;
 245       for (j = 0; j < rows; j++) {
 246         dist[srindex] = galois_single_multiply(tmp, dist[srindex], w);
 247         srindex += cols;
 248       }
 249     }
 250
 251     /* Now, for each element in row i that is not in column 1, you need
 252        to make it zero.  Suppose that this is column j, and the element
 253        at i,j = e.  Then you want to replace all of column j with
 254        (col-j + col-i*e).   Note, that in row i, col-i = 1 and col-j = e.
 255        So (e + 1e) = 0, which is indeed what we want. */
 256
 257     for (j = 0; j < cols; j++) {
 258       tmp = dist[sindex+j];
 259       if (j != i && tmp != 0) {
 260         srindex = j;
 261         siindex = i;
 262         for (k = 0; k < rows; k++) {
 263           dist[srindex] = dist[srindex] ^ galois_single_multiply(tmp, dist[siindex], w);
 264           srindex += cols;
 265           siindex += cols;
 266         }
 267       }
 268     }
 269   }
 270   /* We desire to have row k be all ones.  To do that, multiply
 271      the entire column j by 1/dist[k,j].  Then row j by 1/dist[j,j]. */
 272
 273   sindex = cols*cols;
 274   for (j = 0; j < cols; j++) {
 275     tmp = dist[sindex];
 276     if (tmp != 1) {
 277       tmp = galois_single_divide(1, tmp, w);
 278       srindex = sindex;
 279       for (i = cols; i < rows; i++) {
 280         dist[srindex] = galois_single_multiply(tmp, dist[srindex], w);
 281         srindex += cols;
 282       }
 283     }
 284     sindex++;
 285   }
 286
 287   /* Finally, we'd like the first column of each row to be all ones.  To
 288      do that, we multiply the row by the inverse of the first element. */
 289
 290   sindex = cols*(cols+1);
 291   for (i = cols+1; i < rows; i++) {
 292     tmp = dist[sindex];
 293     if (tmp != 1) {
 294       tmp = galois_single_divide(1, tmp, w);
 295       for (j = 0; j < cols; j++) dist[sindex+j] = galois_single_multiply(dist[sindex+j], tmp, w);
 296     }
 297     sindex += cols;
 298   }
 299
 300   return dist;
 301 }
 302