1 /* -*-c-*-

2 *

3 * $Id: mpmul.h,v 1.2 2004/04/08 01:36:15 mdw Exp $

4 *

5 * Multiply many small numbers together

6 *

7 * (c) 2000 Straylight/Edgeware

8 */

10 /*----- Licensing notice --------------------------------------------------*

11 *

12 * This file is part of Catacomb.

13 *

14 * Catacomb is free software; you can redistribute it and/or modify

15 * it under the terms of the GNU Library General Public License as

16 * published by the Free Software Foundation; either version 2 of the

17 * License, or (at your option) any later version.

18 *

19 * Catacomb is distributed in the hope that it will be useful,

20 * but WITHOUT ANY WARRANTY; without even the implied warranty of

21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the

22 * GNU Library General Public License for more details.

23 *

24 * You should have received a copy of the GNU Library General Public

25 * License along with Catacomb; if not, write to the Free

26 * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,

27 * MA 02111-1307, USA.

28 */

30 #ifndef CATACOMB_MPMUL_H

31 #define CATACOMB_MPMUL_H

33 #ifdef __cplusplus

35 #endif

37 /*----- Header files ------------------------------------------------------*/

39 #ifndef CATACOMB_MP_H

41 #endif

43 /*----- Magic numbers -----------------------------------------------------*/

45 /* --- How the algorithm works --- *

46 *

47 * Multiplication on large integers is least wasteful when the numbers

48 * multiplied are approximately the same size. When a new multiplier is

49 * added to the system, we push it onto a stack. Then we `reduce' the stack:

50 * while the value on the top of the stack is not shorter than the value

51 * below it, replace the top two elements by their product.

52 *

53 * Let %$b$% be the radix of our multiprecision integers, and let %$Z$% be

54 * the maximum number of digits. Then the largest integer we can represent

55 * is %$M - 1 = b^Z - 1$%. We could assume that all of the integers we're

56 * given are about the same size. This would give us the same upper bound as

57 * that derived in `mptext.c'.

58 *

59 * However, we're in less control over our inputs. In particular, if a

60 * sequence of integers with strictly decreasing lengths is input then we're

61 * sunk. Suppose that the stack contains, from top to bottom, %$b^i$%,

62 * %$b^{i+1}$%, ..., %$b^n$%. The final product will therefore be

63 * %$p = b^{(n+i)(n-i+1)/2}$%. We must now find the maximum stack depth

64 * %$d = n - i$% such that %$p > M$%.

65 *

66 * Taking logs of both sides gives that %$(d + 2 i)(d + 1) > 2 Z$%. We can

67 * maximize %$d$% by taking %$i = 0$%, which gives that %$d^2 + d > 2 Z$%, so

68 * %$d$% must be approximately %$(\sqrt{8 Z + 1} - 1)/2$%, which is

69 * uncomfortably large.

70 *

71 * We compromise by choosing double the `mptext' bound and imposing high- and

72 * low-water marks for forced reduction.

73 */

75 #define MPMUL_DEPTH (2 * (CHAR_BIT * sizeof(size_t) + 10))

77 #define HWM (MPMUL_DEPTH - 20)

78 #define LWM (MPMUL_DEPTH / 2)

80 /*----- Data structures ---------------------------------------------------*/

87 #define MPMUL_INIT { 0 }

89 /*----- Functions provided ------------------------------------------------*/

91 /* --- @mpmul_init@ --- *

92 *

93 * Arguments: @mpmul *b@ = pointer to multiplier context to initialize

94 *

95 * Returns: ---

96 *

97 * Use: Initializes a big multiplier context for use.

98 */

102 /* --- @mpmul_add@ --- *

103 *

104 * Arguments: @mpmul *b@ = pointer to multiplier context

105 * @mp *x@ = the next factor to multiply in

106 *

107 * Returns: ---

108 *

109 * Use: Contributes another factor to the mix. It's important that

110 * the integer lasts at least as long as the multiplication

111 * context; this sort of rules out @mp_build@ integers.

112 */

116 /* --- @mpmul_done@ --- *

117 *

118 * Arguments: @mpmul *b@ = pointer to big multiplication context

119 *

120 * Returns: The product of all the numbers contributed.

121 *

122 * Use: Returns a (large) product of numbers. The context is

123 * deallocated.

124 */

128 /* --- @mp_factorial@ --- *

129 *

130 * Arguments: @unsigned long i@ = number whose factorial should be

131 * computed.

132 *

133 * Returns: The requested factorial.

134 */

138 /*----- That's all, folks -------------------------------------------------*/

140 #ifdef __cplusplus

141 }

142 #endif

144 #endif