X-Git-Url: https://git.distorted.org.uk/~mdw/catacomb/blobdiff_plain/ba6e6b64033b1f9de49feccb5c9cd438354481f7..0f00dc4c8eb47e67bc0f148c2dd109f73a451e0a:/gf-arith.c diff --git a/gf-arith.c b/gf-arith.c deleted file mode 100644 index 5a7b3f29..00000000 --- a/gf-arith.c +++ /dev/null @@ -1,309 +0,0 @@ -/* -*-c-*- - * - * $Id$ - * - * Basic arithmetic on binary polynomials - * - * (c) 2004 Straylight/Edgeware - */ - -/*----- Licensing notice --------------------------------------------------* - * - * This file is part of Catacomb. - * - * Catacomb is free software; you can redistribute it and/or modify - * it under the terms of the GNU Library General Public License as - * published by the Free Software Foundation; either version 2 of the - * License, or (at your option) any later version. - * - * Catacomb is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU Library General Public License for more details. - * - * You should have received a copy of the GNU Library General Public - * License along with Catacomb; if not, write to the Free - * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, - * MA 02111-1307, USA. - */ - -/*----- Header files ------------------------------------------------------*/ - -#include "gf.h" - -/*----- Macros ------------------------------------------------------------*/ - -#define MAX(x, y) ((x) >= (y) ? (x) : (y)) - -/*----- Main code ---------------------------------------------------------*/ - -/* --- @gf_add@ --- * - * - * Arguments: @mp *d@ = destination - * @mp *a, *b@ = sources - * - * Returns: Result, @a@ added to @b@. - */ - -mp *gf_add(mp *d, mp *a, mp *b) -{ - MP_DEST(d, MAX(MP_LEN(a), MP_LEN(b)), (a->f | b->f) & MP_BURN); - gfx_add(d->v, d->vl, a->v, a->vl, b->v, b->vl); - d->f = (a->f | b->f) & MP_BURN; - MP_SHRINK(d); - return (d); -} - -/* --- @gf_mul@ --- * - * - * Arguments: @mp *d@ = destination - * @mp *a, *b@ = sources - * - * Returns: Result, @a@ multiplied by @b@. - */ - -mp *gf_mul(mp *d, mp *a, mp *b) -{ - a = MP_COPY(a); - b = MP_COPY(b); - - if (MP_LEN(a) <= MPK_THRESH || MP_LEN(b) <= GFK_THRESH) { - MP_DEST(d, MP_LEN(a) + MP_LEN(b), a->f | b->f | MP_UNDEF); - gfx_mul(d->v, d->vl, a->v, a->vl, b->v, b->vl); - } else { - size_t m = MAX(MP_LEN(a), MP_LEN(b)); - mpw *s; - MP_DEST(d, 2 * m, a->f | b->f | MP_UNDEF); - s = mpalloc(d->a, 3 * m); - gfx_kmul(d->v, d->vl, a->v, a->vl, b->v, b->vl, s, s + 3 * m); - mpfree(d->a, s); - } - - d->f = (a->f | b->f) & MP_BURN; - MP_SHRINK(d); - MP_DROP(a); - MP_DROP(b); - return (d); -} - -/* --- @gf_sqr@ --- * - * - * Arguments: @mp *d@ = destination - * @mp *a@ = source - * - * Returns: Result, @a@ squared. - */ - -mp *gf_sqr(mp *d, mp *a) -{ - MP_COPY(a); - MP_DEST(d, 2 * MP_LEN(a), a->f & MP_BURN); - gfx_sqr(d->v, d->vl, a->v, a->vl); - d->f = a->f & MP_BURN; - MP_SHRINK(d); - MP_DROP(a); - return (d); -} - -/* --- @gf_div@ --- * - * - * Arguments: @mp **qq, **rr@ = destination, quotient and remainder - * @mp *a, *b@ = sources - * - * Use: Calculates the quotient and remainder when @a@ is divided by - * @b@. The destinations @*qq@ and @*rr@ must be distinct. - * Either of @qq@ or @rr@ may be null to indicate that the - * result is irrelevant. (Discarding both results is silly.) - * There is a performance advantage if @a == *rr@. - */ - -void gf_div(mp **qq, mp **rr, mp *a, mp *b) - { - mp *r = rr ? *rr : MP_NEW; - mp *q = qq ? *qq : MP_NEW; - - /* --- Set the remainder up right --- */ - - b = MP_COPY(b); - a = MP_COPY(a); - if (r) - MP_DROP(r); - r = a; - MP_DEST(r, MP_LEN(b) + 2, a->f | b->f); - - /* --- Fix up the quotient too --- */ - - r = MP_COPY(r); - MP_DEST(q, MP_LEN(r), r->f | MP_UNDEF); - MP_DROP(r); - - /* --- Perform the calculation --- */ - - gfx_div(q->v, q->vl, r->v, r->vl, b->v, b->vl); - - /* --- Sort out the sign of the results --- * - * - * If the signs of the arguments differ, and the remainder is nonzero, I - * must add one to the absolute value of the quotient and subtract the - * remainder from @b@. - */ - - q->f = (r->f | b->f) & MP_BURN; - r->f = (r->f | b->f) & MP_BURN; - - /* --- Store the return values --- */ - - MP_DROP(b); - - if (!qq) - MP_DROP(q); - else { - MP_SHRINK(q); - *qq = q; - } - - if (!rr) - MP_DROP(r); - else { - MP_SHRINK(r); - *rr = r; - } -} - -/* --- @gf_irreduciblep@ --- * - * - * Arguments: @mp *f@ = a polynomial - * - * Returns: Nonzero if the polynomial is irreducible; otherwise zero. - */ - -int gf_irreduciblep(mp *f) -{ - unsigned long m; - mp *u = MP_TWO; - mp *v = MP_NEW; - - if (MP_ZEROP(f)) - return (0); - else if (MP_LEN(f) == 1) { - if (f->v[0] < 2) return (0); - if (f->v[0] < 4) return (1); - } - m = (mp_bits(f) - 1)/2; - while (m) { - u = gf_sqr(u, u); - gf_div(0, &u, u, f); - v = gf_add(v, u, MP_TWO); - gf_gcd(&v, 0, 0, v, f); - if (!MP_EQ(v, MP_ONE)) break; - m--; - } - MP_DROP(u); - MP_DROP(v); - return (!m); -} - -/*----- Test rig ----------------------------------------------------------*/ - -#ifdef TEST_RIG - -static int verify(const char *op, mp *expect, mp *result, mp *a, mp *b) -{ - if (!MP_EQ(expect, result)) { - fprintf(stderr, "\n*** %s failed", op); - fputs("\n*** a = ", stderr); mp_writefile(a, stderr, 16); - fputs("\n*** b = ", stderr); mp_writefile(b, stderr, 16); - fputs("\n*** result = ", stderr); mp_writefile(result, stderr, 16); - fputs("\n*** expect = ", stderr); mp_writefile(expect, stderr, 16); - fputc('\n', stderr); - return (0); - } - return (1); -} - -#define RIG(name, op) \ - static int t##name(dstr *v) \ - { \ - mp *a = *(mp **)v[0].buf; \ - mp *b = *(mp **)v[1].buf; \ - mp *r = *(mp **)v[2].buf; \ - mp *c = op(MP_NEW, a, b); \ - int ok = verify(#name, r, c, a, b); \ - mp_drop(a); mp_drop(b); mp_drop(c); mp_drop(r); \ - assert(mparena_count(MPARENA_GLOBAL) == 0); \ - return (ok); \ - } - -RIG(add, gf_add) -RIG(mul, gf_mul) -RIG(exp, gf_exp) - -#undef RIG - -static int tsqr(dstr *v) -{ - mp *a = *(mp **)v[0].buf; - mp *r = *(mp **)v[1].buf; - mp *c = MP_NEW; - int ok = 1; - c = gf_sqr(MP_NEW, a); - ok &= verify("sqr", r, c, a, MP_ZERO); - mp_drop(a); mp_drop(r); mp_drop(c); - assert(mparena_count(MPARENA_GLOBAL) == 0); - return (ok); -} - -static int tdiv(dstr *v) -{ - mp *a = *(mp **)v[0].buf; - mp *b = *(mp **)v[1].buf; - mp *q = *(mp **)v[2].buf; - mp *r = *(mp **)v[3].buf; - mp *c = MP_NEW, *d = MP_NEW; - int ok = 1; - gf_div(&c, &d, a, b); - ok &= verify("div(quotient)", q, c, a, b); - ok &= verify("div(remainder)", r, d, a, b); - mp_drop(a); mp_drop(b); mp_drop(c); mp_drop(d); mp_drop(r); mp_drop(q); - assert(mparena_count(MPARENA_GLOBAL) == 0); - return (ok); -} - -static int tirred(dstr *v) -{ - mp *a = *(mp **)v[0].buf; - int r = *(int *)v[1].buf; - int c = gf_irreduciblep(a); - int ok = 1; - if (r != c) { - ok = 0; - fprintf(stderr, "\n*** irred failed"); - fputs("\n*** a = ", stderr); mp_writefile(a, stderr, 16); - fprintf(stderr, "\n*** r = %d\n", r); - fprintf(stderr, "*** c = %d\n", c); - } - mp_drop(a); - assert(mparena_count(MPARENA_GLOBAL) == 0); - return (ok); -} - -static test_chunk tests[] = { - { "add", tadd, { &type_mp, &type_mp, &type_mp, 0 } }, - { "mul", tmul, { &type_mp, &type_mp, &type_mp, 0 } }, - { "sqr", tsqr, { &type_mp, &type_mp, 0 } }, - { "div", tdiv, { &type_mp, &type_mp, &type_mp, &type_mp, 0 } }, - { "exp", texp, { &type_mp, &type_mp, &type_mp, 0 } }, - { "irred", tirred, { &type_mp, &type_int, 0 } }, - { 0, 0, { 0 } }, -}; - -int main(int argc, char *argv[]) -{ - sub_init(); - test_run(argc, argv, tests, SRCDIR "/tests/gf"); - return (0); -} - -#endif - -/*----- That's all, folks -------------------------------------------------*/