Double factorial and primorial

bodrato at mail.dm.unipi.it bodrato at mail.dm.unipi.it
Mon Dec 19 18:47:47 CET 2011


Ciao,

Il Lun, 19 Dicembre 2011 10:14 am, Niels Möller ha scritto:
> bodrato at mail.dm.unipi.it writes:

> Since the established name is "double factorial", I think one should
> use some reasonable abbreviation of that term. _dblfac_ seems good
> enough to me, if _double_fac_ is too long.

Finally I tested times for _dblfac_:
 - time needed to compute (2n)!! equals time needed for n!, as expected, -
computing (2n+1)!! require 50% more time.

> I don't know how to compute
> log(promordial(n)), but I imagine it will also be pretty large.

I read from http://en.wikipedia.org/wiki/Chebyshev_function :
"the primorial x# is asymptotically equal to exp((1+o(1))x)"

On shell.gmplib.org I tested n=805306366
 time(n!)   log2(n!)   time(n#)  log2(n#)  time(n!!)  log2(n!!)  t(n+1!!)
2.5*10^12  2.2*10^10  2.5*10^11  1.1*10^9  1.1*10^12  1.1*10^10  1.7*10^12

> But for computing the *number* of primes up to, output size is no issue.

Right.

> Hence, the comparison with factoring, 64 bit numbers seem easy.
> Counting the number of primes up to 2^64 may well be a lot harder.
> Maybe there's no better way to compute the number of primes than to
> use the sieve of Erathostenes' to generate them all?

An asymptotically faster sieve than Erathosteses' exists, but I do not
know if it is worth implementing it; as it is O(N/loglogN) compared to
O(N)...

To have the exact count, I believe the only way is generating all primes
(except, of course, tabulating values for 2^i, then generating only the
primes in the shorter interval: ]2^i,n] or ]n,2^(i+1)[.)

In fac_ui.c there already is a function for sieving: bitwise_primesieve.
It fills an mp_ptr "with the characteristic function of composite"s, and
it returns the count of "prime integers in the interval [4, n]". This
means that the following function can count primes in [1, n] (with n>4).

unsigned long primes_up_to_n (unsigned long n)
{
  mp_limb_t *sieve;
  unsigned long result;
  TMP_DECL;

  ASSERT (n > 4);
  TMP_MARK;
  sieve = TMP_ALLOC_LIMBS (primesieve_size (n));
  result = bitwise_primesieve (sieve, n);
  TMP_FREE;

  return result + 2;
}

Timing on shell.gmplib.org with n=2^34:
n           time        result
17179869184 1.3*10^11   762939111

The time is measured in cycles, you can compare this value with the value
for primorial... counting the primes in [1,2^34] required half the time of
(sieving and) multiplying all primes in [1,3*2^28].

I hope the result is correct :-)

Regards,
m

-- 
http://bodrato.it/software/combinatorics.html






More information about the gmp-devel mailing list