/* sieve.c (emx+gcc) */

/* Define BITS to use bit encoding instead of byte encoding */

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#define FALSE   0
#define TRUE    1

#define N(x) ((x)/3)

typedef unsigned long ULONG;
typedef unsigned char BYTE;

#if defined (BITS)

#define index(c) ((c)>>3)
#define mask(c)  (1<<((c)&0x07))
#define BIT_ON(map,c) (map[index(c)] & mask(c))
#define BIT_OFF(map,c) (!BIT_ON(map,c))
#define SET_TRUE(map,c) map[index(c)] |= mask(c)
#define SET_FALSE(map,c) map[index(c)] &= ~mask(c)
#define BYTES(s) (((s)+7)>>3)

#else

#define BIT_ON(map,c) map[c]
#define BIT_OFF(map,c) !map[c]
#define SET_TRUE(map,c) map[c] = 1
#define SET_FALSE(map,c) map[c] = 0
#define BYTES(s) (s)

#endif

static ULONG isqrt (ULONG x);
static void usage (void);


static ULONG isqrt (ULONG x)
{
  ULONG l, r, m;

  l = 1; r = x;
  while (l < r)
    {
      m = (l+r) / 2;
      if (m > 46340) m = 46340;  /* avoid overflow when computing m*m */
      if (m*m < x)
        l = m+1;
      else if (m*m > x)
        r = m-1;
      else
        return m;
    }
  return (l+r)/2;
}


static void usage (void)
{
  fputs ("Usage: sieve [-p] <number>\n", stderr);
  exit (1);
}


int main (int argc, char *argv[])
{
  ULONG sqrt_size, i, j, j_max, n, primes, size;
  int i_incr, j_incr, idx, print_flag;
  size_t bytes;
  char *p;
  BYTE *bitmap;

  print_flag = FALSE;
  for (idx = 1; idx < argc; ++idx)
    if (strcmp (argv[idx], "-p") == 0)
      print_flag = TRUE;
    else
      break;
  if (argc - idx != 1)
    usage ();
  errno = 0;
  size = strtoul (argv[idx], &p, 10);
  if (errno != 0 || *p != 0 || size < 10 || size > 2000000000)
    usage ();
  bytes = BYTES (N (size) + 1);
  bitmap = malloc (bytes);
  if (bitmap == NULL)
    {
      fputs ("Out of memory\n", stderr);
      exit (2);
    }
  memset (bitmap, 0, bytes);
  sqrt_size = isqrt (size);

  n = 1; i_incr = 4;
  for (i = 5; i <= sqrt_size; i += i_incr)
    {
      if (BIT_OFF (bitmap, n))
        {
          j_incr = i_incr; j_max = size / i;
          for (j = i; j <= j_max; j += j_incr)
            {
              SET_TRUE (bitmap, N (i * j));
              j_incr = 6 - j_incr;
            }
        }
      ++n; i_incr = 6 - i_incr;
    }

  if (print_flag)
    printf ("2 3 ");
  n = 1; primes = 2; i_incr = 4;
  for (i = 5; i <= size; i += i_incr)
    {
      if (BIT_OFF (bitmap, n))
        {
          ++primes;
          if (print_flag)
            printf ("%lu ", i);
        }
      ++n; i_incr = 6 - i_incr;
    }
  printf ("\n%lu primes\n", primes);
  return 0;
}
