1// (C) 1996-1998 Luigi Rizzo (luigi@iet.unipi.it) 2// 2009-2010 Jack Lloyd (lloyd@randombit.net) 3// 2011 Billy Brumley (billy.brumley@aalto.fi) 4// 2016-2017 Vivint, Inc. (jeff.wendling@vivint.com) 5// 6// Portions derived from code by Phil Karn (karn@ka9q.ampr.org), 7// Robert Morelos-Zaragoza (robert@spectra.eng.hawaii.edu) and Hari 8// Thirumoorthy (harit@spectra.eng.hawaii.edu), Aug 1995 9// 10// Redistribution and use in source and binary forms, with or without 11// modification, are permitted provided that the following conditions are 12// met: 13// 14// 1. Redistributions of source code must retain the above copyright 15// notice, this list of conditions and the following disclaimer. 16// 17// 2. Redistributions in binary form must reproduce the above copyright 18// notice, this list of conditions and the following disclaimer in the 19// documentation and/or other materials provided with the 20// distribution. 21// 22// THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND ANY EXPRESS OR 23// IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED 24// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE 25// DISCLAIMED. IN NO EVENT SHALL THE AUTHORS BE LIABLE FOR ANY DIRECT, 26// INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES 27// (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR 28// SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 29// HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, 30// STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING 31// IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE 32// POSSIBILITY OF SUCH DAMAGE. 33 34package infectious 35 36import ( 37 "bytes" 38 "errors" 39) 40 41type pivotSearcher struct { 42 k int 43 ipiv []bool 44} 45 46func newPivotSearcher(k int) *pivotSearcher { 47 return &pivotSearcher{ 48 k: k, 49 ipiv: make([]bool, k), 50 } 51} 52 53func (p *pivotSearcher) search(col int, matrix []byte) (int, int, error) { 54 if p.ipiv[col] == false && matrix[col*p.k+col] != 0 { 55 p.ipiv[col] = true 56 return col, col, nil 57 } 58 59 for row := 0; row < p.k; row++ { 60 if p.ipiv[row] { 61 continue 62 } 63 64 for i := 0; i < p.k; i++ { 65 if p.ipiv[i] == false && matrix[row*p.k+i] != 0 { 66 p.ipiv[i] = true 67 return row, i, nil 68 } 69 } 70 } 71 72 return 0, 0, errors.New("pivot not found") 73} 74 75func swap(a, b *byte) { 76 tmp := *a 77 *a = *b 78 *b = tmp 79} 80 81// TODO(jeff): matrix is a K*K array, row major. 82func invertMatrix(matrix []byte, k int) error { 83 pivot_searcher := newPivotSearcher(k) 84 indxc := make([]int, k) 85 indxr := make([]int, k) 86 id_row := make([]byte, k) 87 88 for col := 0; col < k; col++ { 89 icol, irow, err := pivot_searcher.search(col, matrix) 90 if err != nil { 91 return err 92 } 93 94 if irow != icol { 95 for i := 0; i < k; i++ { 96 swap(&matrix[irow*k+i], &matrix[icol*k+i]) 97 } 98 } 99 100 indxr[col] = irow 101 indxc[col] = icol 102 pivot_row := matrix[icol*k:][:k] 103 c := pivot_row[icol] 104 105 if c == 0 { 106 return errors.New("singular matrix") 107 } 108 109 if c != 1 { 110 c = gf_inverse[c] 111 pivot_row[icol] = 1 112 mul_c := gf_mul_table[c][:] 113 114 for i := 0; i < k; i++ { 115 pivot_row[i] = mul_c[pivot_row[i]] 116 } 117 } 118 119 id_row[icol] = 1 120 if !bytes.Equal(pivot_row, id_row) { 121 p := matrix 122 for i := 0; i < k; i++ { 123 if i != icol { 124 c = p[icol] 125 p[icol] = 0 126 addmul(p[:k], pivot_row, c) 127 } 128 p = p[k:] 129 } 130 } 131 132 id_row[icol] = 0 133 } 134 135 for i := 0; i < k; i++ { 136 if indxr[i] != indxc[i] { 137 for row := 0; row < k; row++ { 138 swap(&matrix[row*k+indxr[i]], &matrix[row*k+indxc[i]]) 139 } 140 } 141 } 142 return nil 143} 144 145func createInvertedVdm(vdm []byte, k int) { 146 if k == 1 { 147 vdm[0] = 1 148 return 149 } 150 151 b := make([]byte, k) 152 c := make([]byte, k) 153 154 c[k-1] = 0 155 for i := 1; i < k; i++ { 156 mul_p_i := gf_mul_table[gf_exp[i]][:] 157 for j := k - 1 - (i - 1); j < k-1; j++ { 158 c[j] ^= mul_p_i[c[j+1]] 159 } 160 c[k-1] ^= gf_exp[i] 161 } 162 163 for row := 0; row < k; row++ { 164 index := 0 165 if row != 0 { 166 index = int(gf_exp[row]) 167 } 168 mul_p_row := gf_mul_table[index][:] 169 170 t := byte(1) 171 b[k-1] = 1 172 for i := k - 2; i >= 0; i-- { 173 b[i] = c[i+1] ^ mul_p_row[b[i+1]] 174 t = b[i] ^ mul_p_row[t] 175 } 176 177 mul_t_inv := gf_mul_table[gf_inverse[t]][:] 178 for col := 0; col < k; col++ { 179 vdm[col*k+row] = mul_t_inv[b[col]] 180 } 181 } 182} 183