User talk:Anders Hellström/Pari

Next prime

Note that A151800 can be implemented as nextprime(n+1) which is faster and simpler than your current np() program.

Charles R Greathouse IV 16:35, 12 September 2015 (UTC)

Missing functionality

saving program/session to file.

You might be interested in writebin, which (when given no object) dumps all user variables (including functions) to the specified file. You can read this back to a clean session to get back to where you were, aside from the history entries %1, %2, etc. You might also be interested in default(histfile, 1).

vector(size,i=-1,some_expr(i))

You've already noticed that you can do vector(size,i,some_expr(i-2)). You might also like apply(some_expr, [start..end]).

listing all userdefined functions.

Just type ?0.

Charles R Greathouse IV 16:15, 14 September 2015 (UTC)

indexof

An alternative to

indexof(n)=i=1;while(prime(i)<n,i++);i \\Index of smallest prime > n

would be

indexof(n)=primepi(n-1)+1

which is significantly faster for large n. For example,

for(n=1,1e5,indexof(n))

takes 1min, 5s with your code but only 31ms with mine.

Charles R Greathouse IV 18:24, 16 September 2015 (UTC)

Thoughts

Just some thoughts on your functions:

divides(n,m)=my(r=m/n);r==floor(r)

is most efficiently written

divides(n,m)=m%n == 0

and if you have the need to divide two numbers without remainder x\y is often easier/faster than floor(x/y).

I agree that issemiprime would be very welcome. The PARI library has ways to do partial factorization, so that (for example)

moebius(randomprime(10^100)*randomprime(10^100)*9)

returns 0 almost instantly even though it couldn't reasonably factor such a large number. But this isn't exposed in GP so you either need to factor the number fully (like in your code) or you need to cobble together rudimentary partial factorization in GP, which is slow compared to what the library can do. I have PARI code which does this efficiently, plus a GP script which does the best that it reasonably can without access to the PARI internals.

vecprod1(v,n)=my(p=1);for(i=1,n,p*=v[i]);p; \\product up to n-th element of v

could be

vecprod1(v,n)=factorback(v[1..n]); \\product up to n-th element of v

I'm sure your Ackermann function is not designed to be efficient, but if you wanted it to be:

A(m,n)=if(m>3,if(n,A(m-1,A(m,n-1)),A(m-1,1)),m==3,2^(n+3)-3,m==2,2*n+3,m==1,n+2,n+1)

will happily compute A(4,2) > 10^19728 in under a millisecond.

I think that perhaps

expr(v)={my(l = #v,r=0);for(i=1,l,r += prime(i)^v[i];)}

is intending to be

expr(v)=sum(i=1,#v, prime(i)^v[i])

which could be done somewhat more efficiently (not that it's needed here) as

expr(v)=my(s); forprime(p=2,prime(#v), s += p^v[i]); s

I hope that some of these are useful or interesting.

Charles R Greathouse IV 02:37, 12 December 2016 (UTC)

Thanks for your comments // Anders Hellström 03:58, 19 December 2016 (UTC)