A059113 - OEIS

%I #16 Dec 11 2019 10:51:27

%S 31,187131,347529,658749,2106853,3969147,4657471,4980801,6540807,

%T 10185273,20421219,25690413,42827421,44136183,51598911,73854183,

%U 202765017,299457613,307464339,308953281,314858271,504081669,516979281,600877641,613602549

%N Numbers that are the sum of the divisors of two distinct coprime squares.

%H Donovan Johnson, <a href="/A059113/b059113.txt">Table of n, a(n) for n = 1..1995</a> (terms < 10^16)

%e 31 = sigma(4^2) = sigma(5^2); 187131 = sigma(326^2) = sigma(407^2); etc.

%o (PARI) m=10^5; m2=m^2; s=vector(m); ss=vector(m); v=vector(49); for(j=1, m, s[j]=sigma(j^2)); ss=vecsort(s); c=0; for(j=2, m, if(s[j]>m2, next); if(s[j]>=ss[j], for(k=j, m-1, if(s[j]==ss[k], if(s[j]<>ss[k+1], next(2)); k=m-1)), forstep(k=j-1, 2, -1, if(s[j]==ss[k], if(s[j]<>ss[k-1], next(2)); k=2))); for(k=j+1, m, if(s[j]==s[k], if(gcd(j,k)==1, c++; v[c]=s[j])))); v=vecsort(v); n=0; for(j=2, 49, if(v[j-1]<>v[j], n++; print(n " " v[j]))) /* _Donovan Johnson_, Apr 20 2013 */

%Y Cf. A000203.

%K nonn

%O 1,1

%A _David W. Wilson_, Jan 04 2001