まず，基本的なことから．

もし，奇素数 p だけを考えるならば，

  • p = x^2 + y^2 <==> p = 1 (mod 4)
  • p = x^2 + 2 * y^2 <==> p = 1, 3 (mod 8)
  • p = x^2 + 3 * y^2 <==> p = 1 (mod 6)
  • p = x^2 + 7 * y^2 <==> p = 1, 9, 11 (mod 14)

となる．(これは平方剰余から導出できる．)

よって，4つの表し方が可能な奇素数は

p = 1, 25， 121 (mod 168)

である．これだけなら，賢いふるいを使えば，できる．

問題は合成数．条件が複雑(これで，はじめは brute force に向かったわけだが)．

合成数 n を考える．

以下 p_i はそれぞれの剰余の条件を満たす素数，q_iはそれ以外の素数としておく．

  • n = x^2 + y^2 ==> 素因数として，2か p = 1 (mod 4) が必要．また，2^(2k+1) (q_1^m_1…)^2 か 2^2k p_1^n_1.. (q_1^m_1..)^2 のような形である．
  • n = x^2 + 2 * y^2 ==> 素因数として，p = 1, 3 (mod 8) が必要．また， 2^k p_1^n1… (q_1^m_1…)^2 のような形である．
  • n = x^2 + 3 * y^2 ==> 素因数として，2を2個か， p = 1 (mod 6)が必要．また， 2^2k 3^l (q_1^m_1…)^2 か 2^k 3^l p_1^n1… (q_1^m_1…)^2 のような形である．
  • n = x^2 + 7 * y^2 ==> 素因数として，2を3か4個か， p = 1, 9, 11 (mod 14)が必要．また，2^(4+2k) 7^l (q_1^m_1…)^2 か 2^(3+2k) 7^l (q_1^m_1…)^2 か 2^2 7^l p_1^n_1.. (q_1^m_1…)^2 のような形である．

これで，わけわからなくなったわけですよ．

でも，これから，実は求めたい数は

K^2, K^2*p, K^2*p*q, K^2*p*q*r, K^2*p*q*r*s, …

という形をしていて，p,q,r,sなどは p = 1, 25, 121 (mod 168) を満たす素数であることがわかる．

はじめは，これらの数をまとめて求めようとしていたが，K^2とそれ以外に分けて考えれば良いことに気がついた．

K について探索すれば良いので，探索範囲は √N でうれしい．

(K^2*p のタイプは sqrt( N / p )個ある，K^2*p*q 等も同様)

K^2 に関する条件を整理すると…

  • K^2 = x^2 + y^2 → K が素因数として p = 1 (mod 4) を持つ．
  • K^2 = x^2 + 2 * y^2 → K が素因数として p = 1, 3 (mod 8) を持つ．
  • K^2 = x^2 + 3 * y^2 → K が素因数として 2 か p = 1 (mod 6) を持つ．
  • K^2 = x^2 + 7 * y^2 → K が素因数として 2を2個 か p = 1, 9 11 (mod 14) を持つ．

これならチェックは簡単．

素因数分解を行なうのが気になるが，ふるいで √N までの素因数は持っているので大丈夫．

しかも，探索範囲は √N までだし．