1 /* Copyright (c) 1984, 1986, 1987, 1988, 1989 AT&T */ 2 /* All Rights Reserved */ 3 4 5 /* 6 * Copyright (c) 1980 Regents of the University of California. 7 * All rights reserved. The Berkeley software License Agreement 8 * specifies the terms and conditions for redistribution. 9 */ 10 /* Portions Copyright(c) 1988, Sun Microsystems Inc. */ 11 /* All Rights Reserved */ 12 13 /* 14 * Copyright (c) 1997, by Sun Microsystems, Inc. 15 * All rights reserved. 16 */ 17 18 #ident "%Z%%M% %I% %E% SMI" /* SVr4.0 1.1 */ 19 20 /* LINTLIBRARY */ 21 22 #include <mp.h> 23 #include <sys/types.h> 24 #include "libmp.h" 25 26 static void m_mult(MINT *, MINT *, MINT *); 27 28 void 29 mp_mult(MINT *a, MINT *b, MINT *c) 30 { 31 struct mint x, y; 32 int sign; 33 34 _mp_mcan(a); 35 _mp_mcan(b); 36 if (a->len == 0 || b->len == 0) { 37 _mp_xfree(c); 38 return; 39 } 40 sign = 1; 41 x.len = y.len = 0; 42 _mp_move(a, &x); 43 _mp_move(b, &y); 44 if (a->len < 0) { 45 x.len = -x.len; 46 sign = -sign; 47 } 48 if (b->len < 0) { 49 y.len = -y.len; 50 sign = -sign; 51 } 52 _mp_xfree(c); 53 if (x.len < y.len) { 54 m_mult(&x, &y, c); 55 } else { 56 m_mult(&y, &x, c); 57 } 58 if (sign < 0) 59 c->len = -c->len; 60 if (c->len == 0) 61 _mp_xfree(c); 62 _mp_xfree(&x); 63 _mp_xfree(&y); 64 } 65 66 /* 67 * Knuth 4.3.1, Algorithm M 68 */ 69 static void 70 m_mult(MINT *a, MINT *b, MINT *c) 71 { 72 int i, j; 73 int sum; 74 short bcache; 75 short *aptr; 76 short *bptr; 77 short *cptr; 78 short fifteen = 15; 79 int alen; 80 int blen; 81 82 #define BASEBITS (8 * (unsigned int)sizeof (short) - 1) 83 #define BASE (1 << BASEBITS) 84 #define LOWBITS (BASE - 1) 85 86 alen = a->len; 87 blen = b->len; 88 89 c->len = alen + blen; 90 c->val = _mp_xalloc(c->len, "m_mult"); 91 92 aptr = a->val; 93 bptr = b->val; 94 cptr = c->val; 95 96 sum = 0; 97 bcache = *bptr++; 98 for (i = alen; i > 0; i--) { 99 sum += *aptr++ * bcache; 100 *cptr++ = (short)(sum & LOWBITS); 101 if (sum >= 0) 102 sum >>= fifteen; 103 else 104 sum = 0xfffe0000 | (sum >> fifteen); 105 } 106 *cptr = (short)sum; 107 aptr -= alen; 108 cptr -= alen; 109 cptr++; 110 111 for (j = blen - 1; j > 0; j--) { 112 sum = 0; 113 bcache = *bptr++; 114 for (i = alen; i > 0; i--) { 115 sum += *aptr++ * bcache + *cptr; 116 *cptr++ = (short)(sum & LOWBITS); 117 if (sum >= 0) 118 sum >>= fifteen; 119 else 120 sum = 0xfffe0000 | (sum >> fifteen); 121 } 122 *cptr = (short)sum; 123 aptr -= alen; 124 cptr -= alen; 125 cptr++; 126 } 127 if (c->val[c->len-1] == 0) { 128 c->len--; 129 } 130 } 131