A008849 Numbers n such that the sum of divisors of n^3 is a square.

1, 7, 751530, 4730879, 5260710, 33116153, 37200735, 187062910, 226141311, 259109835, 260405145, 370049418, 522409465, 836308083, 1105725765, 1309440370, 1343713507, 1582989177, 1609505430, 1813768845, 2590345926, 3039492538, 3656866255

Offset: 1

Author

N. J. A. Sloane

Keywords

• nonn
• nice

Comments

In 1657 Fermat challenged the world to find such numbers. [Dickson, Vol. 1, p. 54]

If n is a term and n is not divisible by 7, then 7*n is a term. - Don Dechman (dondechman_2000(AT)yahoo.com), Mar 26 2008

References

• A. H. Beiler, Recreations in the Theory of Numbers, Dover, NY, 1964, p. 9.
• L. E. Dickson, History of the Theory of Numbers. Carnegie Institute Public. 256, Washington, DC, Vol. 1, 1919; Vol. 2, 1920; Vol. 3, 1923, see vol. 1, p. 54.
• Julian Havil, Gamma, Exploring Euler's Constant, Princeton University Press, Princeton and Oxford, 2003, page 92.
• I. Kaplansky, The challenges of Fermat, Wallis and Ozanam (and several related challenges): I. Fermat's first challenge, Preprint, 2002.

Links

• Donovan Johnson, Table of n, a(n) for n = 1..40 (terms < 10^11)
• Donovan Johnson, 73688 terms > 10^11
• Eric Weisstein's World of Mathematics, Fermat's Divisor Problem.

Crossrefs

Cf. A046872, A008850, A048948.

Programs

Mathematica

max = 10^11; primes = {5, 7, 11, 13, 17, 19, 23, 31, 41, 43, 47, 83, 89, 191, 193, 239, 307, 443, 463, 499, 557, 701, 743, 1087, 1487, 2309, 3583, 4373, 5087, 5507, 5807, 44179}; subs = Select[ Times @@@ Subsets[primes, 7], # < max &] // Sort; f[e2_, e3_, p_] := If[n = 2^e2*3^e3*p; IntegerQ[ Sqrt[ DivisorSigma[1, n^3]]], Print[{2^e2, 3^e3, p}]; Sow[n]]; r = Reap[ Scan[ ((f[0, 0, #]; f[0, 1, #]; f[0, 3, #]; f[1, 0, #]; f[1, 1, #]; f[1, 3, #]; f[3, 0, #]; f[3, 1, #]; f[3, 3, #])& ), subs]][[2, 1]]; Select[r, # < max &] // Union (* Jean-François Alcover, Sep 07 2012, after Donovan Johnson *)

PARI

is(n)=issquare(sigma(n^3)) \\ Charles R Greathouse IV, Jun 20 2013

Python

from functools import reduce
from operator import mul
from sympy import factorint, integer_nthroot
A008849_list, n = [], 1
while n < 10**7:
    fs = factorint(n)
    if integer_nthroot(reduce(mul,((p**(3*fs[p]+1)-1)//(p-1) for p in fs),1),2)[1]:
        A008849_list.append(n)
    n += 1 # Chai Wah Wu, Apr 05 2021

Extensions

David W. Wilson has supplied terms a(4) = 4730879 and beyond and verified completeness up to a(3) = 751530

I. Kaplansky and Will Jagy have verified that there are no other terms below 3.8*10^9

3656866255 added by Don Dechman, Mar 26 2008