A004613 - OEIS
A004613
Numbers that are divisible only by primes congruent to 1 mod 4.
42
1, 5, 13, 17, 25, 29, 37, 41, 53, 61, 65, 73, 85, 89, 97, 101, 109, 113, 125, 137, 145, 149, 157, 169, 173, 181, 185, 193, 197, 205, 221, 229, 233, 241, 257, 265, 269, 277, 281, 289, 293, 305, 313, 317, 325, 337, 349, 353, 365, 373, 377, 389, 397, 401, 409, 421
COMMENTS
Also gives solutions z to x^2+y^2=z^4 with gcd(x,y,z)=1 and x,y,z positive. - John Sillcox (johnsillcox(AT)hotmail.com), Feb 20 2004
A062327(a(n)) = A000005(a(n))^2. (These are the only numbers that satisfy this equation.) - Benedikt Otten, May 22 2013
Numbers that are positive integer divisors of 1 + 4*x^2 where x is a positive integer. - Michael Somos, Jul 26 2013
Numbers k such that there is a "knight's move" of Euclidean distance sqrt(k) which allows the whole of the 2D lattice to be reached. For example, a knight which travels 4 units in any direction and then 1 unit at right angles to the first direction moves a distance sqrt(17) for each move. This knight can reach every square of an infinite chessboard.
Also 1/7 of the area of the n-th largest octagon with angles 3*Pi/4, along the perimeter of which there are only 8 nodes of the square lattice - at its vertices. - Alexander M. Domashenko, Feb 21 2024
REFERENCES
David A. Cox, "Primes of the Form x^2 + n y^2", Wiley, 1989.
FORMULA
Numbers of the form x^2 + y^2 where x is even, y is odd and gcd(x, y) = 1.
MAPLE
isA004613 := proc(n)
local p;
for p in numtheory[factorset](n) do
if modp(p, 4) <> 1 then
return false;
end if;
end do:
true;
end proc:
for n from 1 to 200 do
if isA004613(n) then
printf("%d, ", n) ;
end if;
# second Maple program:
q:= n-> andmap(i-> irem(i[1], 4)=1, ifactors(n)[2]):
MATHEMATICA
ok[1] = True; ok[n_] := And @@ (Mod[#, 4] == 1 &) /@ FactorInteger[n][[All, 1]]; Select[Range[421], ok] (* Jean-François Alcover, May 05 2011 *)
Select[Range[500], Union[Mod[#, 4]&/@(FactorInteger[#][[All, 1]])]=={1}&] (* Harvey P. Dale, Mar 08 2017 *)
PROG
(PARI) for(n=1, 1000, if(sumdiv(n, d, isprime(d)*if((d-1)%4, 1, 0))==0, print1(n, ", ")))
(PARI) is(n)=n%4==1 && factorback(factor(n)[, 1]%4)==1 \\ Charles R Greathouse IV, Sep 19 2016
(Magma) [n: n in [1..500] | forall{d: d in PrimeDivisors(n) | d mod 4 eq 1}]; // Vincenzo Librandi, Aug 21 2012
(Haskell)
a004613 n = a004613_list !! (n-1)
a004613_list = filter (all (== 1) . map a079260 . a027748_row) [1..]