chiark - git - mdw - catacomb/blob - mptext-len.c

   1 /* -*-c-*-
   2  *
   3  * $Id$
   4  *
   5  * Work out length of a number's string representation
   6  *
   7  * (c) 2002 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 /*----- Header files ------------------------------------------------------*/
  31
  32 #include "mp.h"
  33 #include "mptext.h"
  34
  35 /*----- Main code ---------------------------------------------------------*/
  36
  37 /* --- @mptext_len@ --- *
  38  *
  39  * Arguments:   @mp *x@ = number to work on
  40  *              @int r@ = radix the number will be expressed in
  41  *
  42  * Returns:     The number of digits needed to represent the number in the
  43  *              given base.  This will not include space for a leading sign
  44  *              (use @MP_NEGP@ to check that, or just add one on for luck);
  45  *              neither will it add space for a terminating null.  In general
  46  *              the answer will be an overestimate.
  47  */
  48
  49 size_t mptext_len(mp *x, int r)
  50 {
  51   unsigned long b = mp_bits(x);
  52   int s, ss = 2;
  53   size_t n;
  54   unsigned d = 0;
  55
  56   /* --- Huh? --- *
  57    *
  58    * The number of digits is at most %$\lceil b \log 2/\log r \rceil$%.  We
  59    * produce an underestimate of %$\log_2 r = \log r/\log 2$% and divide by
  60    * that.  How?  By linear interpolation between known points on the curve.
  61    * The known points are precisely the powers of 2, so we can find a pair
  62    * efficiently by doubling up.  The log curve is convex, so linear
  63    * interpolation between points on the curve is always an underestimate.
  64    *
  65    * The integer maths here is a bit weird, so here's how it works.  If
  66    * %$s = 2^d$% is the power of 2 below %$r$% then we want to compute
  67    * %$\lceil b/(d + (r - s)/s) \rceil = \lceil (b s)/(s(d - 1) + r \rceil$%
  68    * which is %$\lfloor (r + s (b + d - 1) - 1)/(r + s(d - 1)) \rfloor$%.
  69    * Gluing the whole computation together like this makes the code hard to
  70    * read, but means that there are fewer possibilities for rounding errors
  71    * and thus we get a tighter bound.
  72    */
  73
  74   /* --- Find the right pair of points --- */
  75
  76   if (r < 0) r = -r;
  77   do {
  78     s = ss;
  79     d++;
  80     if (r == s) {
  81       n = (b + (d - 1))/d;
  82       goto done;
  83     }
  84     ss = s << 1;
  85   } while (ss <= r);
  86
  87   /* --- Do the interpolation --- */
  88
  89   n = (r + s*(b + d - 1) - 1)/(r + s*(d - 1));
  90
  91   /* --- Fixups --- */
  92
  93 done:
  94   if (!n)
  95     n = 1;
  96   return (n);
  97 }
  98
  99 /*----- That's all, folks -------------------------------------------------*/