3 * Low-level arithmetic on binary polynomials
5 * (c) 2000 Straylight/Edgeware
8 /*----- Licensing notice --------------------------------------------------*
10 * This file is part of Catacomb.
12 * Catacomb is free software; you can redistribute it and/or modify
13 * it under the terms of the GNU Library General Public License as
14 * published by the Free Software Foundation; either version 2 of the
15 * License, or (at your option) any later version.
17 * Catacomb is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 * GNU Library General Public License for more details.
22 * You should have received a copy of the GNU Library General Public
23 * License along with Catacomb; if not, write to the Free
24 * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
28 /*----- Header files ------------------------------------------------------*/
35 /*----- Main code ---------------------------------------------------------*/
37 /* --- @gfx_add@ --- *
39 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
40 * @const mpw *av, *avl@ = first addend vector base and limit
41 * @const mpw *bv, *bvl@ = second addend vector base and limit
45 * Use: Adds two %$\gf{2}$% polynomials. This is the same as
49 void gfx_add(mpw *dv, mpw *dvl,
50 const mpw *av, const mpw *avl,
51 const mpw *bv, const mpw *bvl)
56 while (av < avl || bv < bvl) {
60 a = (av < avl) ? *av++ : 0;
61 b = (bv < bvl) ? *bv++ : 0;
68 /* --- @gfx_acc@ --- *
70 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
71 * @const mpw *av, *avl@ = addend vector base and limit
75 * Use: Adds the addend into the destination. This is considerably
76 * faster than the three-address add call.
79 void gfx_acc(mpw *dv, mpw *dvl, const mpw *av, const mpw *avl)
92 /* --- @gfx_accshift@ --- *
94 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
95 * @const mpw *av, *avl@ = addend vector base and limit
96 * @size_t n@ = number of bits to shift
100 * Use: Shifts the argument left by %$n$% places and adds it to the
101 * destination. This is a primitive used by multiplication and
105 void gfx_accshift(mpw *dv, mpw *dvl, const mpw *av, const mpw *avl, size_t n)
107 size_t c = n / MPW_BITS;
111 /* --- Sort out the shift amounts --- */
118 gfx_acc(dv, dvl, av, avl);
123 /* --- Sort out vector lengths --- */
131 /* --- Now do the hard work --- */
135 *dv++ ^= MPW((y << n) | (x >> c));
142 /* --- @gfx_mul@ --- *
144 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
145 * @const mpw *av, *avl@ = first argument vector base and limit
146 * @const mpw *bv, *bvl@ = second argument vector base and limit
150 * Use: Does multiplication of polynomials over %$\gf{2}$%.
153 void gfx_mul(mpw *dv, mpw *dvl, const mpw *av, const mpw *avl,
154 const mpw *bv, const mpw *bvl)
164 MPSCAN_INITX(&sc, av, avl);
167 while (bv < bvl && dv < dvl) {
169 for (v = av, vv = dv++; v < avl && vv < dvl; v++) {
184 /* --- @gfx_div@ --- *
186 * Arguments: @mpw *qv, *qvl@ = quotient vector base and limit
187 * @mpw *rv, *rvl@ = dividend/remainder vector base and limit
188 * @const mpw *dv, *dvl@ = divisor vector base and limit
192 * Use: Performs division on polynomials over %$\gf{2}$%.
195 void gfx_div(mpw *qv, mpw *qvl, mpw *rv, mpw *rvl,
196 const mpw *dv, const mpw *dvl)
198 size_t dlen, rlen, qlen;
206 assert(((void)"division by zero in gfx_div", dv < dvl));
207 MPX_BITS(dbits, dv, dvl);
217 rvm = 1 << (MPW_BITS - 1);
218 n = MPW_BITS - (dbits % MPW_BITS);
228 gfx_accshift(rvd, rvl, dv, dvl, n);
232 rvm = 1 << (MPW_BITS - 1);
250 /*----- Test rig ----------------------------------------------------------*/
254 #include <mLib/alloc.h>
255 #include <mLib/dstr.h>
256 #include <mLib/quis.h>
257 #include <mLib/testrig.h>
259 #define ALLOC(v, vl, sz) do { \
261 mpw *_vv = xmalloc(MPWS(_sz)); \
262 mpw *_vvl = _vv + _sz; \
267 #define LOAD(v, vl, d) do { \
268 const dstr *_d = (d); \
270 ALLOC(_v, _vl, MPW_RQ(_d->len)); \
271 mpx_loadb(_v, _vl, _d->buf, _d->len); \
276 #define MAX(x, y) ((x) > (y) ? (x) : (y))
278 static void dumpmp(const char *msg, const mpw *v, const mpw *vl)
283 fprintf(stderr, " %08lx", (unsigned long)*--vl);
287 static int vadd(dstr *v)
298 ALLOC(d, dl, MAX(al - a, bl - b) + 1);
300 gfx_add(d, dl, a, al, b, bl);
301 if (!mpx_ueq(d, dl, c, cl)) {
302 fprintf(stderr, "\n*** vadd failed\n");
305 dumpmp("expected", c, cl);
306 dumpmp(" result", d, dl);
310 xfree(a); xfree(b); xfree(c); xfree(d);
314 static int vmul(dstr *v)
325 ALLOC(d, dl, (al - a) + (bl - b));
327 gfx_mul(d, dl, a, al, b, bl);
328 if (!mpx_ueq(d, dl, c, cl)) {
329 fprintf(stderr, "\n*** vmul failed\n");
332 dumpmp("expected", c, cl);
333 dumpmp(" result", d, dl);
337 xfree(a); xfree(b); xfree(c); xfree(d);
341 static int vdiv(dstr *v)
350 ALLOC(a, al, MPW_RQ(v[0].len) + 2); mpx_loadb(a, al, v[0].buf, v[0].len);
354 ALLOC(qq, qql, al - a);
356 gfx_div(qq, qql, a, al, b, bl);
357 if (!mpx_ueq(qq, qql, q, ql) ||
358 !mpx_ueq(a, al, r, rl)) {
359 fprintf(stderr, "\n*** vdiv failed\n");
360 dumpmp(" divisor", b, bl);
361 dumpmp("expect r", r, rl);
362 dumpmp("result r", a, al);
363 dumpmp("expect q", q, ql);
364 dumpmp("result q", qq, qql);
368 xfree(a); xfree(b); xfree(r); xfree(q); xfree(qq);
372 static test_chunk defs[] = {
373 { "add", vadd, { &type_hex, &type_hex, &type_hex, 0 } },
374 { "mul", vmul, { &type_hex, &type_hex, &type_hex, 0 } },
375 { "div", vdiv, { &type_hex, &type_hex, &type_hex, &type_hex, 0 } },
379 int main(int argc, char *argv[])
381 test_run(argc, argv, defs, SRCDIR"/t/gfx");
387 /*----- That's all, folks -------------------------------------------------*/