Sorry not to be clear.
I'd come across a Python function to calculate primepi and/or sum of
primes.
I translated it into c, J and Pari GP with the sorts of performance
times indicated.
The method rests on these properties:
N(v,p) = N(v,p-1) - (N(v/p, p-1) - N(p-1,p-1))
S(v,p) = S(v,p-1) - p(S(v/p, p-1) - S(p-1,p-1))
where N(v,m) (S(v,m)) is the number (sum) of integers in [2,v] that
remain after
sieving by all primes smaller than or equal to m .
The implementation is clever in that it doesn't need to explicitly
include all integers
in the sieve, but instead uses a couple of arrays of size ~ root(n).
My Lehmer effort is based on Weisstein's description in Wolfram World,
and is
MUCH more complicated and less J-like than what I chose to call primepi.
FYI, I've just rediscovered my note to self on the apparent error in
Weisstein's
description:
" NB. Line 1 is wrong; the limits on i,j should read 1<=i<j<=a "
(ie _not_ 1<=i<=j<=a )
OK?
Mike
PS - sorry that table of (very approx) timings displayed poorly for you.
On 23/03/2016 13:33, 'Jon Hough' via Chat wrote:
Sorry, I'm confused, what's the difference between primepi and lehmer?
Anyway, it seems your J Lehmer implementation calculates pi(1e9) in ~2 seconds.
That doesn't seem bad. I mean it's usable.
--------------------------------------------
On Tue, 3/22/16, Mike Day <[email protected]> wrote:
Subject: Re: [Jchat] primes
To: [email protected]
Date: Tuesday, March 22, 2016, 7:32 PM
fwiw, I messed around
with primepi(n) and sum-of-primes-to(n) a few months
ago. (Sum of primes can be dealt with in a
similar fashion to
primepi). Lehmer
didn't perform well in my attempts to apply
it. The Wolfram page on
Lehmer
used to have an error in
one of the summations, so take care if using
http://mathworld.wolfram.com/LehmersFormula.html
(I'm sorry but I can't see if it's
still wrong right now. Weisstein didn't
acknowledge an attempt to report the problem,
and I naven't retained a
copy of
my message.)
Timings for primepi(n), on a previous
laptop:
+------+----------+-------+-------+--------+-----+-----+-----+----+------+
|Times |n= |1000000|1e7 |1e8
|1e9 |1e10 |1e11 |1e12|1e13 |
+------+----------+-------+-------+--------+-----+-----+-----+----+------+
|c++~ |primepi |1ms |5ms
|14ms |150ms|1.1s |9.4s |86s |13m05s|
+------+----------+-------+-------+--------+-----+-----+-----+----+------+
|J |primepi |1ms
|6ms |20ms |62ms |3.4s |16.6s|92s |8m15s |
+------+----------+-------+-------+--------+-----+-----+-----+----+------+
|J |Lehmer |4ms
|22ms |190ms |2.2s |8m23s|....
|....|.... |
+------+----------+-------+-------+--------+-----+-----+-----+----+------+
|PariGP|my-primepi|100ms |616ms |4.4s
|34.6s|4m27s|.... |....|.... |
+------+----------+-------+-------+--------+-----+-----+-----+----+------+
|PariGP|built-in
|0.004ms|0.327ms|0.0011ms|104ms|27.8s|.... |....|.... |
+------+----------+-------+-------+--------+-----+-----+-----+----+------+
(Apologies if this doesn't
display well - I've tried!)
On the present laptop, with an amd chip:
_1 p: <.1e8
5761455
timer'{.{:primepi <.1e8'
+---------+-------+
|0.0179977|5761455|
+---------+-------+
_1 p:
<.1e10
|limit error
| _1 p:<.10000000000
timer'{.{:primepi <.1e10'
+-----+---------+
|3.857|455052511|
+-----+---------+
As the
comparisons above demonstrate, my attempt at coding Lehmer
in J
performed poorly. I code c
and c++ very naively, not knowing them well,
so it's not surprising that J outperforms
my c++ effort for larger n.
The first row of
the Pari GP times are for my attempt at translating my
J function into its script. Pari GP's
built-in primpi function has is
very fast
for lower n, but degrades sharply for larger values. It
might be
relevant that I initialise GP with
parisize = 8000000, primelimit =
4270000000
The
downside, as others have said, is that an efficient
primepi
function does
NOT
scan all primes on its way to counting or summing them!
Code available if
interested.
Mike
On 22/03/2016 00:31, 'Jon
Hough' via Chat wrote:
> Since p: has
a limit (why is that, by the way? I have written a pi(x)
function using Miessel's formula in another language and
it does grind to a slowdown for values much larger than 1e7,
so that might be the reason. I'm assuming internally J
uses Miessel or Lehmer's formula.),
>
it might be better to "cheat" and just use a
hardcoded value for pi(1e9).
>
> https://primes.utm.edu/howmany.html
(scroll down to table 1)
>
>
> Also, possibly the
worst idea ever, why not iterate from 2 to 1e9 and do prime
tests and keep a record of the previously found prime to
compare the final digits with the next prime, and keep a
record of the comparisons. Probably hopelessly slow, and
inefficient.
>
>
--------------------------------------------
> On Tue, 3/22/16, Raul Miller <[email protected]>
wrote:
>
> Subject: Re: [Jchat]
primes
> To: "Chat
forum" <[email protected]>
> Date: Tuesday, March 22,
2016, 12:14 AM
>
> Oh, I coded wrong. The
> limit error was from p:, but
I should not have been
> feeding it the value 1e9.
>
> Sometimes I wonder how I get
anything
> done...
>
> Thanks,
>
> --
> Raul
>
>
> On Mon, Mar
> 21, 2016 at 11:03 AM, Roger
Hui <[email protected]>
> wrote:
>
> >
> > it should be
straightforward to do the extra credit
> problem
> >
> > I should
> have checked. The extra
credit problem is for 100,000,000
> primes.
> >
> > _1 p:
> _1+2^31
> > 105097564
> > p: 1e8
> >
> 2038074751
> >
> > So the
> limit error is not in p: .
> >
> >
> >
> >
> >
> >
> > On Sun, Mar 20, 2016
> at 4:44 PM, Raul Miller
<[email protected]>
> > wrote:
> >
> > > Rosettacode has a
task now, for this
> item.
> > >
> > > http://rosettacode.org/wiki/Prime_conspiracy#J
> > >
> > > And,
> hypothetically speaking, it
should be straightforward to do
> the
> > extra
> > >
> credit problem. But that
doesn't work. Given:
> > >
> > > dgpairs=: 2
> (,'->',])&":/\
10 | p:
> > >
combine=: ~.@[ ,.~ ' ',.~
> ":@,.@(+//.)
> > >
> > > We can try this:
> >
> >
> > > /:~
combine&;/|:
> (~.;#/.~)@dgpairs@((+
i.)/)"1
> > >
> (1e6*i.1e3),.1e6+999>i.1e3
> > >
> > > |limit error:
dgpairs
> > >
> > > ... but p:
> doesn't work for values
like 1e9 (1 p: works, but not
> p:
> > > itself).
> >
> >
> > > And, for example,
Roger has
> worked out some
ways of dealing with large
> > > primes -- see
looking at
> > > http://code.jsoftware.com/wiki/Essays/Primality_Tests
> -- but we don't
> > have
> > > anything that's
a drop in
> replacement for
the p: monad.
> > >
> > > So this presents
something of a
> problem -
how would we tackle a problem
> >
> > like this?
> > >
> > > (Please feel free
to change forum to
> programming if you've got
working
> >
> code
> > > rather than just
> ideas...).
> > >
> >
> > Thanks,
> > >
> >
> > --
> > > Raul
> >
> >
> > >
> > >
> > >
> > > On Tue, Mar
> 15, 2016 at 7:02 AM, Cliff
Reiter <[email protected]>
> > > wrote:
> > >
> > > > A look at the
frequencies of
> pairs of
last digits of successive primes:
> > > >
> /:~({.,#)/.~2]\10|p:4+i.1e7
> > >
> >
> > > > 1 1 446808
> > > >
> > > > 1
> 3 756071
> > > >
> >
> > > 1 7 769924
> > > >
> > > > 1 9 526953
> >
> > >
> > > > 3 1 593196
> > > >
> > > > 3
> 3 422302
> > > >
> >
> > > 3 7 714795
> > > >
> > > > 3 9 769915
> >
> > >
> > > > 7 1 639383
> > > >
> > > > 7
> 3 681759
> > > >
> >
> > > 7 7 422289
> > > >
> > > > 7 9 756852
> >
> > >
> > > > 9 1 820369
> > > >
> > > > 9
> 3 640076
> > > >
> >
> > > 9 7 593275
> > > >
> > > > 9 9 446032
> >
> > >
> > > > The 1 follows
1 as
> rare as 9 follows 9,
but rarer is 3 follows 3 as
> > rare
> > > > as 7
> follows 7. 9 1 most popular!
Very curious. Probably should
> move to
> > > > JProg'
Best,
> Cliff
> > > >
> >
> > >
> > > >
> > > >
> > > >
> On 3/14/2016 12:03 PM, R.E.
Boss wrote:
> >
> > >
> > > >>
> > > >>
> > >
> >
https://www.quantamagazine.org/20160313-mathematicians-discover-prime-conspiracy/
> > > >> R.E.
Boss
> > > >>
> ----------------------------------------------------------------------
> > > For
> > >
> >> information about J
forums see http://www.jsoftware.com/forums.htm
> >
> > >>
> > > >
> > > >
> ----------------------------------------------------------------------
> > > > For
information about J forums
> see http://www.jsoftware.com/forums.htm
> > > >
> > >
> ----------------------------------------------------------------------
> > > For information
about J forums see http://www.jsoftware.com/forums.htm
> > >
> >
> ----------------------------------------------------------------------
> > For information about J
forums see http://www.jsoftware.com/forums.htm
> >
> ----------------------------------------------------------------------
> For information about J
forums see http://www.jsoftware.com/forums.htm
>
----------------------------------------------------------------------
> For information about J forums see http://www.jsoftware.com/forums.htm
---
This email has been checked for viruses by
Avast antivirus software.
https://www.avast.com/antivirus
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
---
This email has been checked for viruses by Avast antivirus software.
https://www.avast.com/antivirus
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
