libaaruformat 1.0
Aaru Data Preservation Suite - Format Library
Loading...
Searching...
No Matches
gf256.c
Go to the documentation of this file.
1/*
2 * This file is part of the Aaru Data Preservation Suite.
3 * Copyright (c) 2019-2026 Natalia Portillo.
4 *
5 * This library is free software; you can redistribute it and/or modify
6 * it under the terms of the GNU Lesser General Public License as
7 * published by the Free Software Foundation; either version 2.1 of the
8 * License, or (at your option) any later version.
9 *
10 * This library is distributed in the hope that it will be useful, but
11 * WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 * Lesser General Public License for more details.
14 *
15 * You should have received a copy of the GNU Lesser General Public
16 * License along with this library; if not, see <http://www.gnu.org/licenses/>.
17 */
18
27
28#include <stddef.h>
29#include <stdint.h>
30#include <string.h>
31
32#include "gf256.h"
33#include "aaruformat/simd.h"
34
35/* -------------------------------------------------------------------------
36 * Log / anti-log tables for GF(2^8) with polynomial 0x11D
37 * Generator element: 2
38 * ------------------------------------------------------------------------- */
39
41static uint8_t gf256_log_table[256];
43static uint8_t gf256_exp_table[512];
44
46static int gf256_tables_initialized = 0;
47
51static void gf256_init_tables(void)
52{
53 if(gf256_tables_initialized) return;
54
55 unsigned x = 1;
56 for(int i = 0; i < 255; i++)
57 {
58 gf256_exp_table[i] = (uint8_t)x;
59 gf256_exp_table[i+255] = (uint8_t)x; /* wrap-around for easy mod 255 */
60 gf256_log_table[x] = (uint8_t)i;
61
62 /* Multiply by generator 2 in GF(2^8) */
63 x <<= 1;
64 if(x & 0x100) x ^= 0x11D;
65 }
66 gf256_log_table[0] = 0; /* Convention: log(0) = 0, unused in mul since we short-circuit */
67 gf256_exp_table[510] = gf256_exp_table[0]; /* Complete wrap */
68 gf256_exp_table[511] = gf256_exp_table[1];
69
70 gf256_tables_initialized = 1;
71}
72
73/* -------------------------------------------------------------------------
74 * Scalar operations
75 * ------------------------------------------------------------------------- */
76
77uint8_t gf256_mul(uint8_t a, uint8_t b)
78{
79 if(a == 0 || b == 0) return 0;
80 gf256_init_tables();
81 return gf256_exp_table[gf256_log_table[a] + gf256_log_table[b]];
82}
83
84uint8_t gf256_div(uint8_t a, uint8_t b)
85{
86 if(a == 0) return 0;
87 /* b must be non-zero */
88 gf256_init_tables();
89 return gf256_exp_table[gf256_log_table[a] + 255 - gf256_log_table[b]];
90}
91
92uint8_t gf256_inv(uint8_t a)
93{
94 /* a must be non-zero */
95 gf256_init_tables();
96 return gf256_exp_table[255 - gf256_log_table[a]];
97}
98
99/* -------------------------------------------------------------------------
100 * SIMD region multiply-accumulate: dst[i] ^= GF_mul(src[i], coeff)
101 *
102 * Technique: 4-bit nibble decomposition.
103 * For a given coeff, precompute:
104 * low_tbl[i] = GF_mul(i, coeff) for i = 0..15
105 * hi_tbl[i] = GF_mul(i<<4, coeff) for i = 0..15
106 * Then for each byte b:
107 * GF_mul(b, coeff) = low_tbl[b & 0x0F] ^ hi_tbl[b >> 4]
108 * This maps to SIMD shuffle (pshufb / vpshufb / vqtbl1q_u8).
109 * ------------------------------------------------------------------------- */
110
114static void gf256_build_mul_tables(uint8_t coeff, uint8_t low_tbl[16], uint8_t hi_tbl[16])
115{
116 gf256_init_tables();
117 for(int i = 0; i < 16; i++)
118 {
119 low_tbl[i] = gf256_mul((uint8_t)i, coeff);
120 hi_tbl[i] = gf256_mul((uint8_t)(i << 4), coeff);
121 }
122}
123
124/* ---------- Scalar fallback ---------- */
125
126static void gf256_mul_region_scalar(uint8_t *dst, const uint8_t *src, uint8_t coeff, size_t len)
127{
128 uint8_t low_tbl[16], hi_tbl[16];
129 gf256_build_mul_tables(coeff, low_tbl, hi_tbl);
130
131 for(size_t i = 0; i < len; i++)
132 dst[i] ^= low_tbl[src[i] & 0x0F] ^ hi_tbl[src[i] >> 4];
133}
134
135static void gf256_xor_region_scalar(uint8_t *dst, const uint8_t *src, size_t len)
136{
137 size_t i = 0;
138
139 /* Process 8 bytes at a time */
140 for(; i + 8 <= len; i += 8)
141 {
142 uint64_t d, s;
143 memcpy(&d, dst + i, 8);
144 memcpy(&s, src + i, 8);
145 d ^= s;
146 memcpy(dst + i, &d, 8);
147 }
148
149 for(; i < len; i++)
150 dst[i] ^= src[i];
151}
152
153/* ---------- x86 SSSE3 ---------- */
154
155#if defined(__x86_64__) || defined(__amd64) || defined(_M_AMD64) || defined(_M_X64) || \
156 defined(__I386__) || defined(__i386__) || defined(__THW_INTEL) || defined(_M_IX86)
157
158#include <tmmintrin.h> /* SSSE3: _mm_shuffle_epi8 */
159
160SSSE3 static void gf256_mul_region_ssse3(uint8_t *dst, const uint8_t *src, uint8_t coeff, size_t len)
161{
162 uint8_t low_tbl[16], hi_tbl[16];
163 gf256_build_mul_tables(coeff, low_tbl, hi_tbl);
164
165 const __m128i low_v = _mm_loadu_si128((const __m128i *)low_tbl);
166 const __m128i hi_v = _mm_loadu_si128((const __m128i *)hi_tbl);
167 const __m128i mask = _mm_set1_epi8(0x0F);
168
169 size_t i = 0;
170 for(; i + 16 <= len; i += 16)
171 {
172 __m128i s = _mm_loadu_si128((const __m128i *)(src + i));
173 __m128i d = _mm_loadu_si128((const __m128i *)(dst + i));
174 __m128i s_lo = _mm_and_si128(s, mask);
175 __m128i s_hi = _mm_and_si128(_mm_srli_epi16(s, 4), mask);
176 __m128i lo = _mm_shuffle_epi8(low_v, s_lo);
177 __m128i hi = _mm_shuffle_epi8(hi_v, s_hi);
178 __m128i r = _mm_xor_si128(_mm_xor_si128(lo, hi), d);
179 _mm_storeu_si128((__m128i *)(dst + i), r);
180 }
181
182 /* Tail */
183 for(; i < len; i++)
184 dst[i] ^= low_tbl[src[i] & 0x0F] ^ hi_tbl[src[i] >> 4];
185}
186
187SSSE3 static void gf256_xor_region_ssse3(uint8_t *dst, const uint8_t *src, size_t len)
188{
189 size_t i = 0;
190 for(; i + 16 <= len; i += 16)
191 {
192 __m128i d = _mm_loadu_si128((const __m128i *)(dst + i));
193 __m128i s = _mm_loadu_si128((const __m128i *)(src + i));
194 _mm_storeu_si128((__m128i *)(dst + i), _mm_xor_si128(d, s));
195 }
196 for(; i < len; i++) dst[i] ^= src[i];
197}
198
199/* ---------- x86 AVX2 ---------- */
200
201#include <immintrin.h> /* AVX2: _mm256_shuffle_epi8 */
202
203AVX2 static void gf256_mul_region_avx2(uint8_t *dst, const uint8_t *src, uint8_t coeff, size_t len)
204{
205 uint8_t low_tbl[16], hi_tbl[16];
206 gf256_build_mul_tables(coeff, low_tbl, hi_tbl);
207
208 /* Broadcast 16-byte tables to both 128-bit lanes of 256-bit register */
209 const __m128i low_128 = _mm_loadu_si128((const __m128i *)low_tbl);
210 const __m128i hi_128 = _mm_loadu_si128((const __m128i *)hi_tbl);
211 const __m256i low_v = _mm256_broadcastsi128_si256(low_128);
212 const __m256i hi_v = _mm256_broadcastsi128_si256(hi_128);
213 const __m256i mask = _mm256_set1_epi8(0x0F);
214
215 size_t i = 0;
216 for(; i + 32 <= len; i += 32)
217 {
218 __m256i s = _mm256_loadu_si256((const __m256i *)(src + i));
219 __m256i d = _mm256_loadu_si256((const __m256i *)(dst + i));
220 __m256i s_lo = _mm256_and_si256(s, mask);
221 __m256i s_hi = _mm256_and_si256(_mm256_srli_epi16(s, 4), mask);
222 __m256i lo = _mm256_shuffle_epi8(low_v, s_lo);
223 __m256i hi = _mm256_shuffle_epi8(hi_v, s_hi);
224 __m256i r = _mm256_xor_si256(_mm256_xor_si256(lo, hi), d);
225 _mm256_storeu_si256((__m256i *)(dst + i), r);
226 }
227
228 /* Tail: SSSE3 for remaining 16-byte chunks, then scalar */
229 const __m128i low_v2 = low_128;
230 const __m128i hi_v2 = hi_128;
231 const __m128i mask2 = _mm_set1_epi8(0x0F);
232 for(; i + 16 <= len; i += 16)
233 {
234 __m128i s = _mm_loadu_si128((const __m128i *)(src + i));
235 __m128i d = _mm_loadu_si128((const __m128i *)(dst + i));
236 __m128i s_lo = _mm_and_si128(s, mask2);
237 __m128i s_hi = _mm_and_si128(_mm_srli_epi16(s, 4), mask2);
238 __m128i lo = _mm_shuffle_epi8(low_v2, s_lo);
239 __m128i hi = _mm_shuffle_epi8(hi_v2, s_hi);
240 __m128i r = _mm_xor_si128(_mm_xor_si128(lo, hi), d);
241 _mm_storeu_si128((__m128i *)(dst + i), r);
242 }
243
244 for(; i < len; i++)
245 dst[i] ^= low_tbl[src[i] & 0x0F] ^ hi_tbl[src[i] >> 4];
246}
247
248AVX2 static void gf256_xor_region_avx2(uint8_t *dst, const uint8_t *src, size_t len)
249{
250 size_t i = 0;
251 for(; i + 32 <= len; i += 32)
252 {
253 __m256i d = _mm256_loadu_si256((const __m256i *)(dst + i));
254 __m256i s = _mm256_loadu_si256((const __m256i *)(src + i));
255 _mm256_storeu_si256((__m256i *)(dst + i), _mm256_xor_si256(d, s));
256 }
257 for(; i + 16 <= len; i += 16)
258 {
259 __m128i d = _mm_loadu_si128((const __m128i *)(dst + i));
260 __m128i s = _mm_loadu_si128((const __m128i *)(src + i));
261 _mm_storeu_si128((__m128i *)(dst + i), _mm_xor_si128(d, s));
262 }
263 for(; i < len; i++) dst[i] ^= src[i];
264}
265
266/* Forward declarations for CPUID-based detection (from simd.c) */
267int have_ssse3(void);
268int have_avx2(void);
269
270#endif /* x86 */
271
272/* ---------- ARM NEON ---------- */
273
274#if defined(__aarch64__) || defined(_M_ARM64) || defined(__arm__) || defined(_M_ARM)
275
276#include <arm_neon.h>
277
278TARGET_WITH_SIMD static void gf256_mul_region_neon(uint8_t *dst, const uint8_t *src, uint8_t coeff, size_t len)
279{
280 uint8_t low_tbl[16], hi_tbl[16];
281 gf256_build_mul_tables(coeff, low_tbl, hi_tbl);
282
283 const uint8x16_t low_v = vld1q_u8(low_tbl);
284 const uint8x16_t hi_v = vld1q_u8(hi_tbl);
285 const uint8x16_t mask = vdupq_n_u8(0x0F);
286
287 size_t i = 0;
288 for(; i + 16 <= len; i += 16)
289 {
290 uint8x16_t s = vld1q_u8(src + i);
291 uint8x16_t d = vld1q_u8(dst + i);
292 uint8x16_t s_lo = vandq_u8(s, mask);
293 uint8x16_t s_hi = vandq_u8(vshrq_n_u8(s, 4), mask);
294 uint8x16_t lo = vqtbl1q_u8(low_v, s_lo);
295 uint8x16_t hi = vqtbl1q_u8(hi_v, s_hi);
296 uint8x16_t r = veorq_u8(veorq_u8(lo, hi), d);
297 vst1q_u8(dst + i, r);
298 }
299
300 for(; i < len; i++)
301 dst[i] ^= low_tbl[src[i] & 0x0F] ^ hi_tbl[src[i] >> 4];
302}
303
304TARGET_WITH_SIMD static void gf256_xor_region_neon(uint8_t *dst, const uint8_t *src, size_t len)
305{
306 size_t i = 0;
307 for(; i + 16 <= len; i += 16)
308 {
309 uint8x16_t d = vld1q_u8(dst + i);
310 uint8x16_t s = vld1q_u8(src + i);
311 vst1q_u8(dst + i, veorq_u8(d, s));
312 }
313 for(; i < len; i++) dst[i] ^= src[i];
314}
315
316int have_neon(void);
317
318#endif /* ARM */
319
320/* -------------------------------------------------------------------------
321 * Dispatch functions
322 * ------------------------------------------------------------------------- */
323
324void gf256_mul_region(uint8_t *dst, const uint8_t *src, uint8_t coeff, size_t len)
325{
326 if(coeff == 0) return;
327 if(coeff == 1) { gf256_xor_region(dst, src, len); return; }
328
329 gf256_init_tables();
330
331#if defined(__x86_64__) || defined(__amd64) || defined(_M_AMD64) || defined(_M_X64) || \
332 defined(__I386__) || defined(__i386__) || defined(__THW_INTEL) || defined(_M_IX86)
333 if(have_avx2()) { gf256_mul_region_avx2(dst, src, coeff, len); return; }
334 if(have_ssse3()) { gf256_mul_region_ssse3(dst, src, coeff, len); return; }
335#endif
336
337#if defined(__aarch64__) || defined(_M_ARM64) || defined(__arm__) || defined(_M_ARM)
338 if(have_neon()) { gf256_mul_region_neon(dst, src, coeff, len); return; }
339#endif
340
341 gf256_mul_region_scalar(dst, src, coeff, len);
342}
343
344void gf256_xor_region(uint8_t *dst, const uint8_t *src, size_t len)
345{
346#if defined(__x86_64__) || defined(__amd64) || defined(_M_AMD64) || defined(_M_X64) || \
347 defined(__I386__) || defined(__i386__) || defined(__THW_INTEL) || defined(_M_IX86)
348 if(have_avx2()) { gf256_xor_region_avx2(dst, src, len); return; }
349 if(have_ssse3()) { gf256_xor_region_ssse3(dst, src, len); return; }
350#endif
351
352#if defined(__aarch64__) || defined(_M_ARM64) || defined(__arm__) || defined(_M_ARM)
353 if(have_neon()) { gf256_xor_region_neon(dst, src, len); return; }
354#endif
355
356 gf256_xor_region_scalar(dst, src, len);
357}
gf256_div
uint8_t gf256_div(uint8_t a, uint8_t b)
Divide two elements in GF(2^8).
Definition gf256.c:84
gf256_xor_region
void gf256_xor_region(uint8_t *dst, const uint8_t *src, size_t len)
XOR a region: dst[i] ^= src[i] for all i.
Definition gf256.c:344
gf256_inv
uint8_t gf256_inv(uint8_t a)
Compute multiplicative inverse in GF(2^8).
Definition gf256.c:92
gf256_log_table
static uint8_t gf256_log_table[256]
Log table: gf256_log[x] = discrete log base 2 of x in GF(2^8).
Definition gf256.c:41
gf256_tables_initialized
static int gf256_tables_initialized
Flag to ensure tables are initialized exactly once.
Definition gf256.c:46
gf256_init_tables
static void gf256_init_tables(void)
Initialize log/antilog tables for GF(2^8) with polynomial 0x11D.
Definition gf256.c:51
gf256_mul_region_scalar
static void gf256_mul_region_scalar(uint8_t *dst, const uint8_t *src, uint8_t coeff, size_t len)
Definition gf256.c:126
gf256_xor_region_scalar
static void gf256_xor_region_scalar(uint8_t *dst, const uint8_t *src, size_t len)
Definition gf256.c:135
gf256_exp_table
static uint8_t gf256_exp_table[512]
Anti-log (exp) table: gf256_exp[i] = 2^i mod P.
Definition gf256.c:43
gf256_mul_region
void gf256_mul_region(uint8_t *dst, const uint8_t *src, uint8_t coeff, size_t len)
Multiply-accumulate a region: dst[i] ^= GF_mul(src[i], coeff) for all i.
Definition gf256.c:324
gf256_build_mul_tables
static void gf256_build_mul_tables(uint8_t coeff, uint8_t low_tbl[16], uint8_t hi_tbl[16])
Build the two 16-byte nibble lookup tables for a given coefficient.
Definition gf256.c:114
gf256_mul
uint8_t gf256_mul(uint8_t a, uint8_t b)
Multiply two elements in GF(2^8) with polynomial 0x11D.
Definition gf256.c:77