Pierpont primes: Difference between revisions

Line 1,630:

I also tried using a priority queue, popping the smallest (and any duplicates of it) and pushing *2 and *3 on every iteration, which worked pretty well but ended up about 20% slower, with about 1500 untested candidates left in the queue by the time it found the 2000th second kind.

<lang Phix>-- ~~demo/rosetta/Pierpont_primes.exw~~

-- demo/rosetta/Pierpont_primes.exw

include mpfr.e

with javascript_semantics

include mpfr.e

function pierpont(integer n)

sequence p = {{mpz_init(2)},{}},

function pierpont(integer n)

s = {mpz_init(1)}

sequence p = {{mpz_init(2)},{}},

integer i2 = 1, i3 = 1

s = {mpz_init(1)}

mpz {n2, n3, t} = mpz_inits(3)

integer i2 = 1, i3 = 1

randstate state = gmp_randinit_mt()

mpz {n2, n3, t} = mpz_inits(3)

atom t1 = time()+1

atom t1 = time()+1

while min(length(p[1]),length(p[2])) < n do

while min(length(p[1]),length(p[2])) < n do

mpz_mul_si(n2,s[i2],2)

mpz_mul_si(n2,s[i2],2)

mpz_mul_si(n3,s[i3],3)

mpz_mul_si(n3,s[i3],3)

if mpz_cmp(n2,n3)<0 then

if mpz_cmp(n2,n3)<0 then

mpz_set(t,n2)

mpz_set(t,n2)

i2 += 1

i2 += 1

else

~~mpz_set(t,n3)~~

else

mpz_set(t,n3)

i3 += 1

i3 += 1

end if

end if

if mpz_cmp(t,s[$]) > 0 then

if mpz_cmp(t,s[$]) > 0 then

s = append(s, mpz_init_set(t))

s = append(s, mpz_init_set(t))

mpz_add_ui(t,t,1)

mpz_add_ui(t,t,1)

if length(p[1]) < n and mpz_probable_prime_p(t,state) then

if length(p[1]) < n and mpz_prime(t) then

p[1] = append(p[1],mpz_init_set(t))

p[1] = append(p[1],mpz_init_set(t))

end if

end if

mpz_sub_ui(t,t,2)

mpz_sub_ui(t,t,2)

if length(p[2]) < n and mpz_probable_prime_p(t,state) then

if length(p[2]) < n and mpz_prime(t) then

p[2] = append(p[2],mpz_init_set(t))

p[2] = append(p[2],mpz_init_set(t))

end if

end if

end if

end if

if time()>t1 then

if time()>t1 then

printf(1,"Found: 1st:%d/%d, 2nd:%d/%d\r",

printf(1,"Found: 1st:%d/%d, 2nd:%d/%d\r",

{length(p[1]),n,length(p[2]),n})

{length(p[1]),n,length(p[2]),n})

t1 = time()+1

t1 = time()+1

end if

end if

end while

end while

return p

return p

end function

end function

constant limit = 2000 -- 2 mins

--constant limit = ~~1000~~ -- ~~8.1s~~

--constant limit = 2000 -- 2 mins

--constant limit = ~~250~~ -- 0.1s

--constant limit = 1000 -- 8.1s

constant limit = 250 -- 0.1s

atom t0 = time()

atom t0 = time()

sequence p = pierpont(limit)

sequence p = pierpont(limit)

constant fs = {"first","second"}

constant fs = {"first","second"}

for i=1 to length(fs) do

for i=1 to length(fs) do

printf(1,"First 50 Pierpont primes of the %s kind:\n",{fs[i]})

printf(1,"First 50 Pierpont primes of the %s kind:\n",{fs[i]})

for j=1 to 50 do

for j=1 to 50 do

mpfr_printf(1,"%8Zd ", p[i][j])

printf(1,"%8s ", {mpz_get_str(p[i][j])})

if mod(j,10)=0 then printf(1,"\n") end if

if mod(j,10)=0 then printf(1,"\n") end if

end for

end for

printf(1,"\n")

printf(1,"\n")

end for

end for

constant t = {250,1000,2000}

constant t = {250,1000,2000}

for i=1 to length(t) do

for i=1 to length(t) do

integer ti = t[i]

integer ti = t[i]

if ti>limit then exit end if

if ti>limit then exit end if

for j=1 to length(fs) do

for j=1 to length(fs) do

string zs = shorten(mpz_get_str(p[j][ti]))

string zs = shorten(mpz_get_str(p[j][ti]))

printf(1,"%dth Pierpont prime of the %s kind: %s\n",{ti,fs[j],zs})

printf(1,"%dth Pierpont prime of the %s kind: %s\n",{ti,fs[j],zs})

end for

end for

printf(1,"\n")

printf(1,"\n")

end for

end for

printf(1,"Took %s\n", elapsed(time()-t0))</lang>

printf(1,"Took %s\n", elapsed(time()-t0))

<pre>

Line 1,724:

Line 1,726:

Took 2 minutes and 01s

</pre>

Note that shorten() has recently been added as a builtin, in honour of this and the several other dozen rc

tasks that previously all re-implemented some variation of it. Should you not yet have it, just use this:

<lang Phix>function shorten(string s, what="digits", integer ml=20)

integer l = length(s)

string ls = sprintf(" (%d %s)",{l,what})

if l>ml*2+3+length(ls) then

s[ml..-ml] = "..."

s &= ls

end if

return s

end function</lang>

=={{header|Prolog}}==