|
%I #15 Aug 22 2018 01:54:27
%S 1,2,4,8,11,13,16,19,22,23,26,29,32,37,38,41,43,44,46,47,52,53,58,59,
%T 61,64,67,71,74,76,79,82,83,86,88,89,92,94,97,101,103,104,106,109,113,
%U 116,118,121,122,128,131,134,137,139,142,143,148,149,151,152,157
%N Numbers having no divisor d > 1 that is a binary palindrome (i.e., an element of A006995).
%e The nonunit divisors of 22 are 2,11,22 and none of these are binary palindromes.
%p dmax:= 10: # to get all terms with at most dmax binary digits
%p N:= 2^dmax-1:
%p revdigs:= proc(n)
%p local L, Ln, i;
%p L:= convert(n, base, 2);
%p Ln:= nops(L);
%p add(L[i]*2^(Ln-i), i=1..Ln);
%p end proc:
%p P:= {}:
%p for d from 2 to dmax do
%p if d::even then
%p P:= P union {seq(2^(d/2)*x + revdigs(x), x=2^(d/2-1)..2^(d/2)-1)}
%p else
%p m:= (d-1)/2;
%p B:={seq(2^(m+1)*x + revdigs(x), x=2^(m-1)..2^m-1)};
%p P:= P union B union map(`+`, B, 2^m)
%p fi
%p od:
%p L:= Vector(N,1):
%p for t in P do
%p L[[seq(k,k=t..N,t)]]:= 0
%p od:
%p select(t -> L[t]=1, [$1..N]); # _Robert Israel_, Aug 21 2018
%o (PARI) isok(n) = #select(x->((binary(x) == Vecrev(binary(x))) && (x>1)), divisors(n)) == 0; \\ _Michel Marcus_, Aug 21 2018
%Y Cf. A006995.
%K nonn,base
%O 1,2
%A _Jeffrey Shallit_, Aug 21 2018
|
