%I #18 Jun 22 2024 15:48:37

%S 2,6,14,18,38,42,50,54,62,74,86,114,122,126,134,146,150,158,162,186,

%T 206,218,222,254,258,266,302,314,326,342,350,366,378,386,398,402,422,

%U 434,438,450,458,474,482,518,542,554,558,566,578,602,618,626,654,662,666,674,686,734,746,758,762,774,794

%N Integers k such that sigma(sigma(2*k))=2*sigma(sigma(k)); sigma=A000203.

%H Graeme L. Cohen and H. J. J. te Riele, <a href="https://doi.org/10.1080/10586458.1996.10504580">Iterating the sum-of-divisors function</a>, Exp. Math., 5 (1996), 91-100.

%p with(numtheory):

%p P := proc (q)

%p local n, result:

%p result := []:

%p for n to q do

%p if sigma(sigma(2*n)) = 2*sigma(sigma(n)) then

%p result := [op(result), n]:

%p end if

%p end do:

%p print(result):

%p end proc:

%p P(10^3);

%t Select[Range[800],DivisorSigma[1,DivisorSigma[1,2#]]==2DivisorSigma[1,DivisorSigma[1,#]]&] (* _Stefano Spezia_, Jun 05 2024 *)

%o (Python)

%o from sympy import divisor_sigma as sigma

%o def P(q):

%o result = []

%o for n in range(1, q + 1):

%o if sigma(sigma(2 * n)) == 2 * sigma(sigma(n)):

%o result.append(n)

%o print(result)

%o P(10**3)

%Y Cf. A000203, A051027.

%K nonn

%O 1,1

%A _Rafik Khalfi_, Jun 04 2024