bug#12350: Composites identified as primes in factor.c (when HAVE_GMP)

[Torbjorn Granlund]

  >   #endif
  >   #define LONGLONG_STANDALONE     /* Don't require GMP's longlong.h mdep 
files */
  >   #define ASSERT(x)               /* FIXME make longlong.h really 
standalone */
  > + #define __GMP_DECLSPEC          /* FIXME make longlong.h really 
standalone */
  
  Oh, that must be used in an #else branch that I'm not exercising.
  It triggers a "symbol defined but never used" warning.
  Thanks.
  
__GMP_DECLSPEC is used from longlong.h when count_leading_zeros ask for
__clz_tab.

  >   #define __clz_tab factor_clz_tab /* Rename to avoid glibc collision */
  
  BTW, why does __clz_tab use the "__" prefix in the first place?
  
It makes sense in glibc.  In GMP we prepend something like __gmp_.  Here
in a non-library we prepend factor_.  I'd like to keep the __ for src
compatibility.

  >   # endif
  >   # include "longlong.h"
  >   # ifdef COUNT_LEADING_ZEROS_NEED_CLZ_TAB
  > ! const unsigned char factor_clz_tab[129] =
  >   {
  >     1,2,3,3,4,4,4,4,5,5,5,5,5,5,5,5,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
  >     7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
  > --- 144,150 ----
  >   #endif
  >   #include "longlong.h"
  >   #ifdef COUNT_LEADING_ZEROS_NEED_CLZ_TAB
  > ! static const unsigned char factor_clz_tab[129] =
  
  Thanks.  Done.
  
Actually, that patch is wrong.  :-(

In longlong.h, we declare it without static.  This causes a type clash.

I suggest that we kep clz_tab non-static.  Do you have any other
suggestions?

  >     while (mpz_cmp_ui (n, 1) != 0)
  >       {
  > ***************
  > *** 2287,2293 ****
  >
  >   /* Read white-space delimited items. Return 1 on success, 0 on EOF.
  >      Exit on I/O errors. */
  > ! int
  >   read_item (struct inbuf *bufstruct)
  >   {
  >     int c;
  > --- 2288,2294 ----
  >
  >   /* Read white-space delimited items. Return 1 on success, 0 on EOF.
  >      Exit on I/O errors. */
  > ! static int
  
  I didn't bother with this one, because in coreutils
  I've just switched back to using readtokens to read input.
  
Ok.  Hopefully, that does not cause too much slowdown for large ranges.
With this fast factoring code, input is a significant part of the time!

I have a variant of input the code, which replaces all the input code
with a single function which freads 4 KiB blocks and plays with
sentinels.  Here is its non-GMP part.  I left this code out for a 20%
performance penalty.

void
read_stdin_and_factor ()
{
  uintmax_t n1, n0;
  unsigned int carry;
  unsigned int c, tc;
  char *bp, *be;
  int valid_char_seen;
  size_t nread;
  enum { BUFSIZE = 4096 };
  char buf[BUFSIZE + 1];

  buf[0] = '\377';                      /* sentinel */
  bp = buf;
  be = buf;

 restart:
  n0 = 0;
  valid_char_seen = 0;

#ifndef OPTIMISE_SINGLE_WORD
#define OPTIMISE_SINGLE_WORD
#endif

#if OPTIMISE_SINGLE_WORD
  /* Loop while we have one word */
  while (LIKELY (n0 < (~(uintmax_t) 0 - 9) / 10))
    {
      c = (unsigned char) *bp++;
      tc = c - '0';

      if (UNLIKELY (tc >= 10))
        {
          if (bp > be)          /* did we hit the sentinel? */
            {
              nread = fread (buf, 1, BUFSIZE, stdin);
              if (nread == 0)
                {
                  if (valid_char_seen)
                    {
                      factor_and_report (0, n0);
                    }
                  return;
                }

              buf[nread] = '\377';      /* sentinel */
              bp = buf;
              be = buf + nread;
              continue;
            }
          else
            {
              if (valid_char_seen)
                {
                  factor_and_report (0, n0);
                  valid_char_seen = 0;
                }
              if (! isspace (c))
                {
                  fprintf (stderr, "`%c' is not a valid positive integer\n", c);
                }
            }
          goto restart;
        }

      /* got valid digit */
      n0 = 10 * n0 + tc;
      valid_char_seen = 1;
    }
#endif

  n1 = 0;

  /* Loop while we have two words */
  for (;;)
    {
      c = (unsigned char) *bp++;
      tc = c - '0';

      if (UNLIKELY (tc >= 10))
        {
          if (bp > be)          /* did we hit the sentinel? */
            {
              nread = fread (buf, 1, BUFSIZE, stdin);
              if (nread == 0)
                {
                  if (valid_char_seen)
                    {
                      factor_and_report (n1, n0);
                    }
                  return;
                }

              buf[nread] = '\377';      /* sentinel */
              bp = buf;
              be = buf + nread;
              continue;
            }
          else
            {
              if (valid_char_seen)
                {
                  factor_and_report (n1, n0);
                  valid_char_seen = 0;
                }
              if (! isspace (c))
                {
                  fprintf (stderr, "`%c' is not a valid positive integer\n", c);
                }
            }
          goto restart;
        }

      /* got valid digit */
      carry = (n0 >> (W_TYPE_SIZE - 3)) + (n0 >> (W_TYPE_SIZE - 1));
      carry += 10 * n0 < 2 * n0;

      if (UNLIKELY (n1 > ~(uintmax_t)0 / 2 / 10))
        break;          /* overflow, error unless we have GMP */

      n1 = 10 * n1;
      n0 = 10 * n0 + tc;
      carry += n0 < tc;
      n1 += carry;

      if (UNLIKELY (((n1 << 1) >> 1) != n1))
        break;          /* overflow, error unless we have GMP */

      valid_char_seen = 1;
    }
#if HAVE_GMP
  // FIXME: keep going, non-trivial, since any buffer might be insufficient
  // * Could balk and do O(n^2) algo instead of calling GMP's inp function
  // * Or we could let the buffer size determine, grab chunks from it, feed
  //   such chunks to GMP functions, multiply by suitable powers.  Good, tricky.
  // * Or skip this nifty function altogether when HAVE_GMP?  Bad!
  abort ();
#else
  fprintf (stderr, "number is too large\n");
  // FIXME: perhaps skip over rest of number
#endif
}


-- 
Torbjörn