Karim Belabas on Wed, 30 May 2012 23:27:22 +0200 |
[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]
forprime() |
Hi pari-dev, adapting an implementation provided by Charles Greathouse (thanks!), I just committed to 'master' a general iterator on arbitrary ranges of primes, which is not limited by 'primelimit'. It uses a variety of strategies: the private primetable up to primelimit, a sieve up to primelimit^2, and nextprime() from then on. \\ f() compute the sum of primes in [a,b] in a silly way (18:22) gp > f(a, b = a+10^5) = my(s = 0); forprime(p = a, b, s += p); s; (18:22) gp > default(primelimit,2) \\ kill primetable (18:22) gp > f(10^20) time = 152 ms. %2 = 211500000000000105447503 (18:22) gp > f(0, 10^9) time = 9,810 ms. %3 = 24739512092254535 It helps *a little* to have precomputed lots of primes, but not much: (18:23) gp > default(primelimit,10^9) \\ cache a huge primetable time = 2,364 ms. (18:23) gp > f(0,10^9) time = 7,464 ms. %4 = 24739512092254535 Total time (initialization + forprime iteration) is about the same, cached primes help by saving about 25% in later iterations (7.4s vs. 9.8s), not a major saving [ in this case, the loop body is almost trivial; in non trivial cases, the cost of prime generation will be masked by the time needed to process the loop ] Current libpari code still hardcodes in many places the old NEXT_PRIME_VIA_DIFF() machinery (using exclusively the private prime table) and is still limited by 'primelimit': if the prime table is too small, the function raises an exception. As I replace these by the new iterators -- forprime_next() in libpari --, I expect that it will become mostly useless to precompute lots of primes; say beyond 10^6, which already provides fast primes up to (10^6)^2 = 10^12. N.B. There's no support for primes in arithmetic progressions yet [ forprimestep() ? forstepprime() ? ], but we already have support for infinite loops: omitting the upper bound as in forprime(p = 2, /* empty */, ...blah...) runs through "all" primes in succession. Cheers, K.B. -- Karim Belabas, IMB (UMR 5251) Tel: (+33) (0)5 40 00 26 17 Universite Bordeaux 1 Fax: (+33) (0)5 40 00 69 50 351, cours de la Liberation http://www.math.u-bordeaux1.fr/~belabas/ F-33405 Talence (France) http://pari.math.u-bordeaux1.fr/ [PARI/GP] `