git.distorted.org.uk Git - u/mdw/catacomb/blob - mpmul.h

   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  */
   9
  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  */
  29
  30 #ifndef CATACOMB_MPMUL_H
  31 #define CATACOMB_MPMUL_H
  32
  33 #ifdef __cplusplus
  34   extern "C" {
  35 #endif
  36
  37 /*----- Header files ------------------------------------------------------*/
  38
  39 #ifndef CATACOMB_MP_H
  40 #  include "mp.h"
  41 #endif
  42
  43 /*----- Magic numbers -----------------------------------------------------*/
  44
  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  */
  74
  75 #define MPMUL_DEPTH (2 * (CHAR_BIT * sizeof(size_t) + 10))
  76
  77 #define HWM (MPMUL_DEPTH - 20)
  78 #define LWM (MPMUL_DEPTH / 2)
  79
  80 /*----- Data structures ---------------------------------------------------*/
  81
  82 typedef struct mpmul {
  83   size_t i;
  84   mp *v[MPMUL_DEPTH];
  85 } mpmul;
  86
  87 #define MPMUL_INIT { 0 }
  88
  89 /*----- Functions provided ------------------------------------------------*/
  90
  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  */
  99
 100 extern void mpmul_init(mpmul */*b*/);
 101
 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  */
 113
 114 extern void mpmul_add(mpmul */*b*/, mp */*x*/);
 115
 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  */
 125
 126 extern mp *mpmul_done(mpmul */*b*/);
 127
 128 /* --- @mp_factorial@ --- *
 129  *
 130  * Arguments:   @unsigned long i@ = number whose factorial should be
 131  *                      computed.
 132  *
 133  * Returns:     The requested factorial.
 134  */
 135
 136 extern mp *mp_factorial(unsigned long /*i*/);
 137
 138 /*----- That's all, folks -------------------------------------------------*/
 139
 140 #ifdef __cplusplus
 141   }
 142 #endif
 143
 144 #endif