Charles Greathouse on Tue, 18 Sep 2012 21:23:16 +0200


[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]

Re: forprime


> As for me not looping over all composits, but looping over integers
> with a given number of divisors will be of great interest.

I brought up an idea a while back of a sieve-based iterator
"forfactored" (the name could use work, perhaps?) which would make
forfactored(f = 1, N, func1(f); func2(f); func3(f))
a faster version of
for(n = 1, N, f=factor(n); func1(f); func2(f); func3(f))

This would certainly do what you want.

> How do you use it? Could you give one or two concrete examples?

Well, the most recent example would be this program to compute terms of A011774:

sp(n)=my(f=factor(n)); n*prod(i=1, #f[, 1], 1-1/f[i, 1]) + prod(i=1,
#f[, 1], (f[i, 1]^(f[i, 2]+1)-1)/(f[i, 1]-1))
p=2; forprime(q=3, 1e6, for(n=p+1, q-1, if(sp(n)%n==0, print1(n", "))); p=q)

which would be written

sp(n)=my(f=factor(n)); n*prod(i=1, #f[, 1], 1-1/f[i, 1]) + prod(i=1,
#f[, 1], (f[i, 1]^(f[i, 2]+1)-1)/(f[i, 1]-1))
forcomposite(n=4, 1e6, if(sp(n)%n==0, print1(n", ")))

But this is really an argument for forfactored, which (with
appropriate support from sigma and eulerphi) would be written

forfactored(n=4, 1e6, if(!isprime(n) && (sigma(n)+eulerphi(n))%n == 0,
print1(factorback(n)", ")))

which is much nicer than either of the above and just as efficient.

A better example would be A064623:

f(n)=my(f=factor(n)[, 1]); #f[, 1]>1 && n%sum(i=1, #f, f[i])==0
list(lim)=my(v=List(), p=2); forprime(q=p+1, nextprime(lim),
for(n=p+1, q-1, if(f(n), for(i=1, #v, if(n%v[i]==0, next(2)));
listput(v, n))); p=q); Vec(v)

which would be written

f(n)=my(f=factor(n)[, 1]); #f[, 1]>1 && n%sum(i=1, #f, f[i])==0
list(lim)=my(v=List()); for(n=4, lim, if(f(n), for(i=1, #v,
if(n%v[i]==0, next(2))); listput(v, n))

But it's hard to think of examples off the top of my head. Very often
it just comes up in the middle of a problem where I find the need and
I just code something.

Charles Greathouse
Analyst/Programmer
Case Western Reserve University

On Tue, Sep 18, 2012 at 12:29 PM, Karim Belabas
<Karim.Belabas@math.u-bordeaux1.fr> wrote:
> * Charles Greathouse [2012-09-18 15:32]:
>> I found forcomposite useful enough that I wrote a GP script to do it,
>> just for convenience's sake.
>
> How do you use it? Could you give one or two concrete examples?
>
> Thanks !
>
>     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]