greedy computation of pseudoperfect numbers?

Maximilian Hasler Maximilian.Hasler at martinique.univ-ag.fr
Thu Apr 10 01:34:19 CEST 2008


To complete the record, I offer my PARI code (easily adapted from that
of 138000) and some more values:


%N A000001 a(n) = least semiprime such that all subsets of
{a(1),...,a(n)} have a different sum.
%S A000001 4,6,9,14,22,57,111,218,445,879,1754,
 3518,7034,14069,28129,56271,112529,225073,450139,900274
%o A000001 (PARI) {s=1;p=0; for( n=1,20,
 until( !bitand( s, s>>p ), until( 2==bigomega(p++),)); s+=s<<p; print1( p","))}

Maximilian
PS: adding the experimentally (and eventually theoretically) founded
(largely sub-optimal)
 acceleration "p=p*5\3-2" at the end of the loop, these values are
calculated in little more than 1 sec, else calculating bigomega() for
all integers up to 900 000 takes about 7 sec...

NOTE: I noticed before sending my own sequence that the definition
"a(n) = least xxx such that all subsets of {a(1),...,a(n)} have a
different sum."
does not exclude to take all a(n) equal to a(1), that's why I changed
the wording
to "the subsets of {...} sum to 2^n different values"
(or you say explicitely a(n)>a(n-1), but this in turn does not work for n=0.)



On 4/9/08, Jonathan Post <jvospost3 at gmail.com> wrote:
> Thank you, Richard Mathar.  I submitted this as follows.
>
>  NEW SEQUENCE FROM Jonathan Vos Post
>
>  %I A000001
>  %S A000001 4, 6, 9, 14, 22, 57, 111, 218, 445, 879, 1754
>  %N A000001 a(n) = least semiprime such that all subsets of
>  {a(1),...,a(n)} have a different sum.
>  %C A000001 Maple code and extension by R. J. Mathar
>  (mathar(AT)strw.leidenuniv.nl)
>  %D A000001 J. Stanley Benkoski and P. Erdos, On weird and
>  pseudoperfect numbers, Math. Comp. 28, no. 126, 617–623, 1974.
>  %p A000001 isA001358 := proc(n)
>
>       if numtheory[bigomega](n) = 2 then
>               true;
>       else
>               false;
>       fi ;
>  end:
>
>  setsum := proc(S)
>       add(i,i=S) ;
>  end:
>
>  a := [4] :
>  while true do
>       for anxt from op(-1,a)+1 do
>               if isA001358(anxt) then
>                       aset := combinat[powerset](convert(a,set) union {anxt} );
>                       sset := {} ;
>                       for s in aset do
>                               sset := sset union { setsum(s) } ;
>                       od:
>                       if nops(sset) = nops(aset) then
>                               a := [op(a),anxt] ;
>                               print(a) ;
>                               break ;
>                       fi ;
>               fi ;
>       od:
>  od:
>
> %Y A000001 Cf. A001358, A006037, A064934.
>  %O A000001 1,1
>  %K A000001 ,easy,more,nonn,
>  %A A000001 Jonathan Vos Post (jvospost3 at gmail.com), Apr 09 2008
>  RH
>  RA 192.20.225.32
>  RU
>  RI
>
>






More information about the SeqFan mailing list