%I #22 Oct 02 2024 16:08:56
%S 1,3,3,6,3,10,5,10,7,10,7,20,7,18,10,18,9,26,11,20,18,26,13,36,11,26,
%T 19,36,15,36,17,34,26,34,18,52,19,42,26,36,21,68,23,52,26,50,25,68,29,
%U 42,34,52,27,74,26,68,42,58,31,72,31,66,50,66,26,100,35,68,50,68
%N Number of solutions to x + y == x^2 + y^2 (mod n) with x <= y.
%H Alois P. Heinz, <a href="/A376646/b376646.txt">Table of n, a(n) for n = 1..4000</a>
%p a:=proc(n)
%p local x,y,count;
%p count:=0:
%p for x from 0 to n-1 do
%p for y from x to n-1 do
%p if (x+y) mod n =(x^2+y^2) mod n then count:=count+1; fi;
%p od:
%p od:
%p count;
%p end:
%p # second Maple program:
%p a:= n-> add(add(`if`(x^2-x+y^2-y mod n=0, 1, 0), x=0..y), y=0..n-1):
%p seq(a(n), n=1..70); # _Alois P. Heinz_, Oct 01 2024
%o (PARI) a(n) = sum(y=0, n-1, sum(x=0, y, (x+y) % n == (x^2+y^2) % n)); \\ _Michel Marcus_, Oct 01 2024
%o (Python)
%o def A376646(n):
%o c = 0
%o for x in range(n):
%o m = x*(1-x)%n
%o c += sum(1 for y in range(x,n) if y*(y-1)%n==m)
%o return c # _Chai Wah Wu_, Oct 02 2024
%K nonn,look
%O 1,2
%A _W. Edwin Clark_, Sep 30 2024