mul.txt (2461B)
1 Multiplication 2 3 This describes the multiplication algorithm used by the MPI library. 4 5 This is basically a standard "schoolbook" algorithm. It is slow -- 6 O(mn) for m = #a, n = #b -- but easy to implement and verify. 7 Basically, we run two nested loops, as illustrated here (R is the 8 radix): 9 10 k = 0 11 for j <- 0 to (#b - 1) 12 for i <- 0 to (#a - 1) 13 w = (a[j] * b[i]) + k + c[i+j] 14 c[i+j] = w mod R 15 k = w div R 16 endfor 17 c[i+j] = k; 18 k = 0; 19 endfor 20 21 It is necessary that 'w' have room for at least two radix R digits. 22 The product of any two digits in radix R is at most: 23 24 (R - 1)(R - 1) = R^2 - 2R + 1 25 26 Since a two-digit radix-R number can hold R^2 - 1 distinct values, 27 this insures that the product will fit into the two-digit register. 28 29 To insure that two digits is enough for w, we must also show that 30 there is room for the carry-in from the previous multiplication, and 31 the current value of the product digit that is being recomputed. 32 Assuming each of these may be as big as R - 1 (and no larger, 33 certainly), two digits will be enough if and only if: 34 35 (R^2 - 2R + 1) + 2(R - 1) <= R^2 - 1 36 37 Solving this equation shows that, indeed, this is the case: 38 39 R^2 - 2R + 1 + 2R - 2 <= R^2 - 1 40 41 R^2 - 1 <= R^2 - 1 42 43 This suggests that a good radix would be one more than the largest 44 value that can be held in half a machine word -- so, for example, as 45 in this implementation, where we used a radix of 65536 on a machine 46 with 4-byte words. Another advantage of a radix of this sort is that 47 binary-level operations are easy on numbers in this representation. 48 49 Here's an example multiplication worked out longhand in radix-10, 50 using the above algorithm: 51 52 a = 999 53 b = x 999 54 ------------- 55 p = 98001 56 57 w = (a[jx] * b[ix]) + kin + c[ix + jx] 58 c[ix+jx] = w % RADIX 59 k = w / RADIX 60 product 61 ix jx a[jx] b[ix] kin w c[i+j] kout 000000 62 0 0 9 9 0 81+0+0 1 8 000001 63 0 1 9 9 8 81+8+0 9 8 000091 64 0 2 9 9 8 81+8+0 9 8 000991 65 8 0 008991 66 1 0 9 9 0 81+0+9 0 9 008901 67 1 1 9 9 9 81+9+9 9 9 008901 68 1 2 9 9 9 81+9+8 8 9 008901 69 9 0 098901 70 2 0 9 9 0 81+0+9 0 9 098001 71 2 1 9 9 9 81+9+8 8 9 098001 72 2 2 9 9 9 81+9+9 9 9 098001 73 74 ------------------------------------------------------------------ 75 This Source Code Form is subject to the terms of the Mozilla Public 76 # License, v. 2.0. If a copy of the MPL was not distributed with this 77 # file, You can obtain one at http://mozilla.org/MPL/2.0/.