%I #23 Feb 27 2023 15:54:45
%S 53,106,211,212,397,403,417,419,422,424,437,441,459,781,794,801,806,
%T 817,833,834,838,839,841,844,848,865,874,882,885,918,979,1481,1549,
%U 1562,1565,1571,1573,1585,1588,1589,1602,1612,1613,1634,1637,1665,1666,1667,1668
%N Numbers k > 1 such that, in base 2, k and k^2 contain the same digits in the same proportion.
%H Alois P. Heinz, <a href="/A061656/b061656.txt">Table of n, a(n) for n = 1..10000</a>
%e 53 = 110101_2 and 53^2 = 101011111001_2.
%p p:= n-> add(x^i, i=convert(n, base, 2)):
%p a:= proc(n) option remember; local k;
%p for k from 1+`if`(n=1, 0, a(n-1))
%p while p(k)*2<>p(k^2) do od; k
%p end:
%p seq(a(n), n=1..50); # _Alois P. Heinz_, May 10 2015
%t b2pQ[n_]:=Module[{bn=IntegerDigits[n,2],b2n=IntegerDigits[n^2,2], cbn0, cb2n0}, cbn0=Count[bn,0];cb2n0=Count[b2n,0];cbn0>0&&cb2n0>0 && Count[ bn,1]/cbn0==Count[b2n,1]/cb2n0]; Select[Range[1700],b2pQ] (* _Harvey P. Dale_, Jan 25 2012 *)
%o (Python)
%o from fractions import Fraction
%o from itertools import count, islice
%o def f(i, j):
%o bi, bj = bin(i)[2:], bin(j)[2:]
%o pi = [Fraction(bi.count(d), len(bi)) for d in "01"]
%o pj = [Fraction(bj.count(d), len(bj)) for d in "01"]
%o return pi == pj
%o def ok(n): return f(n, n**2)
%o print([k for k in range(2, 1700) if ok(k)]) # _Michael S. Branicky_, Feb 27 2023
%Y Cf. A061657, A061658, A061659, A061660, A061661, A061662, A061663, A114258, A061664.
%K base,easy,nonn
%O 1,1
%A _Erich Friedman_, Jun 16 2001
%E Offset changed to 1 by _Alois P. Heinz_, May 10 2015