A391107 - OEIS

%I #26 Dec 06 2025 09:10:56

%S 65,85,130,145,170,185,195,205,221,255,260,265,290,305,325,340,365,

%T 370,377,390,410,425,435,442,445,455,481,485,493,505,510,520,530,533,

%U 545,555,565,580,585,595,610,615,625,629,650,663,680,685,689,697,715,725,730,740,745,754,765,780

%N Integers k such that the equation x^2 + y^2 = 2*k^2 has at least 4 solutions in positive integers x < y.

%C These numbers are "square-sum-rich" and represent potential centers for a 3 X 3 magic square of distinct squares. Equivalent to n = (x^2 + y^2)/2 with at least 4 distinct pairs (x<y).

%C If k = Product_i p_i^e_i, where p_i are distinct primes, then k is a term if and only if Product_{i:p_i mod 4 == 1} 2*e_i+1 >= 9. If k is a term, then k*m is a term as well. A131574 is a subsequence. - _Chai Wah Wu_, Dec 01 2025

%H Chai Wah Wu, <a href="/A391107/b391107.txt">Table of n, a(n) for n = 1..10000</a>

%o (Python)

%o import math

%o def count_pairs(n):

%o     L = 2 * n**2

%o     count = 0

%o     for x in range(1, int(math.isqrt(L)) + 1):

%o         y2 = L - x*x

%o         y = int(math.isqrt(y2))

%o         if y*y == y2 and y > x:

%o             count += 1

%o     return count

%o # Generate list of n ≤ MAX with ≥4 pairs

%o MAX = 1500

%o vals = [n for n in range(1, MAX+1) if count_pairs(n) >= 4]

%o print(vals)

%o (Python)

%o from math import prod

%o from itertools import count, islice

%o from sympy import factorint

%o def A391107_gen(startvalue=1): # generator of terms >= startvalue

%o     for k in count(max(startvalue,1)):

%o         if prod((e<<1)|1 for p, e in factorint(k).items() if p&3==1)>=9:

%o             yield k

%o A391107_list = list(islice(A391107_gen(),58)) # _Chai Wah Wu_, Dec 01 2025

%Y Cf. A084648, A002654, A001481, A131574.

%K nonn

%O 1,1

%A _Md. Rad Sarar Anando_, Nov 29 2025