/* -*-c-*-
*
- * $Id: mpx-kmul.c,v 1.3 1999/12/13 15:35:01 mdw Exp $
+ * $Id: mpx-kmul.c,v 1.6 2000/10/08 12:11:01 mdw Exp $
*
* Karatsuba's multiplication algorithm
*
/*----- Revision history --------------------------------------------------*
*
* $Log: mpx-kmul.c,v $
+ * Revision 1.6 2000/10/08 12:11:01 mdw
+ * Use @mpx_ueq@ instead of @MPX_UCMP@.
+ *
+ * Revision 1.5 2000/07/29 17:04:02 mdw
+ * Remove useless header `mpscan.h'.
+ *
+ * Revision 1.4 2000/06/17 11:42:11 mdw
+ * Moved the Karatsuba macros into a separate file for better sharing.
+ * Fixed some comments.
+ *
* Revision 1.3 1999/12/13 15:35:01 mdw
* Simplify and improve.
*
#include <stdio.h>
#include "mpx.h"
+#include "mpx-kmac.h"
/*----- Tweakables --------------------------------------------------------*/
# define KARATSUBA_CUTOFF 2
#endif
-/*----- Addition macros ---------------------------------------------------*/
-
-#define UADD(dv, av, avl) do { \
- mpw *_dv = (dv); \
- const mpw *_av = (av), *_avl = (avl); \
- mpw _c = 0; \
- \
- while (_av < _avl) { \
- mpw _a, _b; \
- mpd _x; \
- _a = *_av++; \
- _b = *_dv; \
- _x = (mpd)_a + (mpd)_b + _c; \
- *_dv++ = MPW(_x); \
- _c = _x >> MPW_BITS; \
- } \
- while (_c) { \
- mpd _x = (mpd)*_dv + (mpd)_c; \
- *_dv++ = MPW(_x); \
- _c = _x >> MPW_BITS; \
- } \
-} while (0)
-
-#define UADD2(dv, dvl, av, avl, bv, bvl) do { \
- mpw *_dv = (dv), *_dvl = (dvl); \
- const mpw *_av = (av), *_avl = (avl); \
- const mpw *_bv = (bv), *_bvl = (bvl); \
- mpw _c = 0; \
- \
- while (_av < _avl || _bv < _bvl) { \
- mpw _a, _b; \
- mpd _x; \
- _a = (_av < _avl) ? *_av++ : 0; \
- _b = (_bv < _bvl) ? *_bv++ : 0; \
- _x = (mpd)_a + (mpd)_b + _c; \
- *_dv++ = MPW(_x); \
- _c = _x >> MPW_BITS; \
- } \
- *_dv++ = _c; \
- while (_dv < _dvl) \
- *_dv++ = 0; \
-} while (0)
-
-#define USUB(dv, av, avl) do { \
- mpw *_dv = (dv); \
- const mpw *_av = (av), *_avl = (avl); \
- mpw _c = 0; \
- \
- while (_av < _avl) { \
- mpw _a, _b; \
- mpd _x; \
- _a = *_av++; \
- _b = *_dv; \
- _x = (mpd)_b - (mpd)_a - _c; \
- *_dv++ = MPW(_x); \
- if (_x >> MPW_BITS) \
- _c = 1; \
- else \
- _c = 0; \
- } \
- while (_c) { \
- mpd _x = (mpd)*_dv - (mpd)_c; \
- *_dv++ = MPW(_x); \
- if (_x >> MPW_BITS) \
- _c = 1; \
- else \
- _c = 0; \
- } \
-} while (0)
-
/*----- Main code ---------------------------------------------------------*/
/* --- @mpx_kmul@ --- *
/* --- How the algorithm works --- *
*
- * Let %$A = xb + y$% and %$B = ub + v$%. Then, simply by expanding, %$AB
- * = x u b^2 + b(x v + y u) + y v$%. That's not helped any, because I've
- * got four multiplications, each four times easier than the one I started
- * with. However, note that I can rewrite the coefficient of %$b$% as
- * %$xv + yu = (x + y)(u + v) - xu - yv$%. The terms %$xu$% and %$yv$%
+ * Let %$A = xb + y$% and %$B = ub + v$%. Then, simply by expanding,
+ * %$AB = x u b^2 + b(x v + y u) + y v$%. That's not helped any, because
+ * I've got four multiplications, each four times easier than the one I
+ * started with. However, note that I can rewrite the coefficient of %$b$%
+ * as %$xv + yu = (x + y)(u + v) - xu - yv$%. The terms %$xu$% and %$yv$%
* I've already calculated, and that leaves only one more multiplication to
* do. So now I have three multiplications, each four times easier, and
* that's a win.
#include <mLib/alloc.h>
#include <mLib/testrig.h>
-#include "mpscan.h"
-
#define ALLOC(v, vl, sz) do { \
size_t _sz = (sz); \
mpw *_vv = xmalloc(MPWS(_sz)); \
ALLOC(s, sl, 2 * m + 32);
mpx_kmul(d, dl, a, al, b, bl, s, sl);
- if (MPX_UCMP(d, dl, !=, c, cl)) {
+ if (!mpx_ueq(d, dl, c, cl)) {
fprintf(stderr, "\n*** umul failed\n");
dumpmp(" a", a, al);
dumpmp(" b", b, bl);