136 lines
3.7 KiB
C
136 lines
3.7 KiB
C
/* LibTomPoly, Polynomial Basis Math -- Tom St Denis
|
|
*
|
|
* LibTomPoly is a public domain library that provides
|
|
* polynomial basis arithmetic support. It relies on
|
|
* LibTomMath for large integer support.
|
|
*
|
|
* This library is free for all purposes without any
|
|
* express guarantee that it works.
|
|
*
|
|
* Tom St Denis, tomstdenis@iahu.ca, http://poly.libtomcrypt.org
|
|
*/
|
|
#include <tompoly.h>
|
|
|
|
|
|
/* Divides a by b (uses the char of b as the char of the result) only for
|
|
* polynomials over GF(p^k)[x]
|
|
*/
|
|
int pb_div(pb_poly *a, pb_poly *b, pb_poly *c, pb_poly *d)
|
|
{
|
|
int err, x;
|
|
pb_poly p, q, r;
|
|
mp_int tmp, tmp2;
|
|
|
|
/* if ||b|| > ||a|| it's a simple copy over */
|
|
if (b->used > a->used) {
|
|
if (d != NULL) {
|
|
pb_copy(a, d);
|
|
}
|
|
if (c != NULL) {
|
|
pb_zero(c);
|
|
}
|
|
}
|
|
|
|
/* zero divisor */
|
|
if (b->used == 0) {
|
|
return MP_VAL;
|
|
}
|
|
|
|
/* compare chars */
|
|
if (mp_cmp(&(a->characteristic), &(b->characteristic)) != MP_EQ) {
|
|
return MP_VAL;
|
|
}
|
|
|
|
/* we don't support char zero! */
|
|
if (mp_iszero(&(b->characteristic)) == MP_YES) {
|
|
return MP_VAL;
|
|
}
|
|
|
|
/* get a copy of "a" to work with */
|
|
if ((err = pb_init_copy(&p, a)) != MP_OKAY) {
|
|
return err;
|
|
}
|
|
|
|
/* init a temp quotient */
|
|
if ((err = pb_init_size(&q, &(b->characteristic), a->used - b->used + 1)) != MP_OKAY) {
|
|
goto __P;
|
|
}
|
|
|
|
/* init a temp polynomial we can work with */
|
|
if ((err = pb_init_size(&r, &(b->characteristic), a->used)) != MP_OKAY) {
|
|
goto __Q;
|
|
}
|
|
|
|
/* now init an mp_int we can work with */
|
|
if ((err = mp_init(&tmp)) != MP_OKAY) {
|
|
goto __R;
|
|
}
|
|
|
|
/* and the inverse of the leading digit of b(x) */
|
|
if ((err = mp_init(&tmp2)) != MP_OKAY) {
|
|
goto __TMP;
|
|
}
|
|
if ((err = mp_invmod(&(b->terms[b->used - 1]), &(b->characteristic), &tmp2)) != MP_OKAY) {
|
|
goto __TMP2;
|
|
}
|
|
|
|
/* now let's reduce us some digits */
|
|
for (x = a->used - 1; x >= b->used-1; x--) {
|
|
/* skip zero digits */
|
|
if (mp_iszero(&(p.terms[x])) == MP_YES) {
|
|
continue;
|
|
}
|
|
|
|
/* now we have a leading digit of p(x), call it A
|
|
and a leading digit of b(x), call it B
|
|
|
|
Now we want to reduce it so we need A + CB = 0 which turns into
|
|
A + CB = 0
|
|
CB = -A
|
|
C = -A/B
|
|
|
|
So we multiply b(x) by C * x^k [where k = ||p(x)|| - ||b(x)||], add that to p(x) and
|
|
we must reduced one digit
|
|
*/
|
|
|
|
/* multiply 1/B [tmp2] by A */
|
|
if ((err = mp_mulmod(&tmp2, &(p.terms[x]), &(b->characteristic), &tmp)) != MP_OKAY) { goto __TMP2; }
|
|
|
|
/* tmp is now a term of the quotient */
|
|
if ((err = mp_copy(&tmp, &(q.terms[x - b->used + 1]))) != MP_OKAY) { goto __TMP2; }
|
|
|
|
/* create r(x) = C */
|
|
pb_zero(&r);
|
|
if ((err = mp_copy(&tmp, &(r.terms[0]))) != MP_OKAY) { goto __TMP2; }
|
|
r.used = 1;
|
|
|
|
/* now multiply r(x) by b(x)*x^k and subtract from p(x) */
|
|
if ((err = pb_mul(b, &r, &r)) != MP_OKAY) { goto __TMP2; }
|
|
if ((err = pb_lshd(&r, x - b->used + 1)) != MP_OKAY) { goto __TMP2; }
|
|
if ((err = pb_sub(&p, &r, &p)) != MP_OKAY) { goto __TMP2; }
|
|
}
|
|
|
|
/* setup q properly (so q.used is valid) */
|
|
q.used = a->used;
|
|
pb_clamp(&q);
|
|
|
|
/* ok so now p(x) is the remainder and q(x) is the quotient */
|
|
if (c != NULL) {
|
|
pb_exch(&q, c);
|
|
}
|
|
if (d != NULL) {
|
|
pb_exch(&p, d);
|
|
}
|
|
|
|
err = MP_OKAY;
|
|
|
|
__TMP2: mp_clear(&tmp2);
|
|
__TMP : mp_clear(&tmp);
|
|
__R : pb_clear(&r);
|
|
__Q : pb_clear(&q);
|
|
__P : pb_clear(&p);
|
|
return err;
|
|
}
|
|
|
|
|