４ｎ＋１型素数の性質　 　　　　　　　　　　　　　　　　　　　　　　

　次の表は、１０００以下の素数表である。その中で、４ｎ＋１型の素数は赤字の数である。
全部で８０個（１６８個中）ある。

　このページでは、４ｎ＋１型素数が持ついろいろな性質をまとめてみたい。

　素数が無限個あるように、４ｎ＋１型素数も無限個ある。また、（偶数）²＋１を素因数分解

すると、その素因数はすべて４ｎ＋１型素数である。また、逆に、４ｎ＋１型素数は、ある偶数

に対して、（偶数）²＋１の素因数となる。

例　ｘ²＋１ の形の整数の素因数は、２ または ４ｎ＋１型の素数であることを確認せよ。

　まず、いくつか計算してみよう。

　ｘ＝１　のとき、　ｘ²＋１＝２

　ｘ＝２　のとき、　ｘ²＋１＝５　→　５ は、４ｎ＋１型の素数

　ｘ＝３　のとき、　ｘ²＋１＝１０＝２・５　→　５ は、４ｎ＋１型の素数

　ｘ＝４　のとき、　ｘ²＋１＝１７　→　１７ は、４ｎ＋１型の素数

　ｘ＝５　のとき、　ｘ²＋１＝２６＝２・１３　→　１３ は、４ｎ＋１型の素数

　ｘ＝６　のとき、　ｘ²＋１＝３７　→　３７ は、４ｎ＋１型の素数

　ｘ＝７　のとき、　ｘ²＋１＝５０＝２・５・５　→　５ は、４ｎ＋１型の素数

　ｘ＝８　のとき、　ｘ²＋１＝６５＝５・１３　→　５、１３ はともに、４ｎ＋１型の素数

　ｘ＝９　のとき、　ｘ²＋１＝８２＝２・４１　→　４１ は、４ｎ＋１型の素数

　ｘ＝１０　のとき、　ｘ²＋１＝１０１　→　１０１ は、４ｎ＋１型の素数

　　　　・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・

　上記の計算で気がつくことは、

　ｘ が偶数のとき、ｘ²＋１ の形の整数の素因数は、４ｎ＋１型の素数のみ

になっていることだろう。

　ｘ が偶数のとき、ｘ²＋１ は、４で割って１余る数で、素因数 ２ は持たず奇数である。

　ｘ が奇数のとき、ｘ²＋１ は、素因数 ２ を必ず持つ。

　そこで、

　　数 ｘ²＋１ の奇数の素因数は、４ｎ＋１型の素数のみ

であることは、次のように証明される。

（証明）　ｘ＝１、２　のとき、命題は成立している。

　３以上の整数 ｘ に対して、初めて命題が成立しないものとする。

　ｘ が偶数のとき

　　ｘ²＋１＝Ａ・Ｂ　（Ａ は、４ｎ＋３型の素数で、Ｂ は整数）と分解されるものとする。

　このとき、左辺は、４で割って１余る数なので、Ｂは、４で割って３余る数でなければならな

い。よって、適当に、Ｂの素因数のうち ４ｎ＋３型の素数とＡとを入れ替えることにより、

Ａ≦Ｂ　としても一般性を失わない。

　ここで、　ｘ＜Ａ　とすると、　ｘ＋１≦Ａ≦Ｂ　で、

　　　ｘ²＋１＝Ａ・Ｂ≧（ｘ＋１）²＝ｘ²＋２ｘ＋１

となり、矛盾するので、　Ａ≦ｘ　である。

　このとき、　（ｘ−Ａ）²＋１＝ｘ²＋１−２ｘＡ＋Ａ²＝Ａ・Ｂ−２ｘＡ＋Ａ²＝Ａ（Ａ＋Ｂ−２ｘ）

　ｘ−Ａ＜ｘ　で、帰納法の仮定より、（ｘ−Ａ）²＋１ の奇数の素因数は、４ｎ＋１型の素数の

みのはず。ところが、奇数 Ａ は、４ｎ＋３型の素数なので矛盾。

　したがって、３以上の偶数 ｘ に対して、命題は成立する。

　次に、ｘ が奇数のとき

　　ｘ²＋１＝２・Ｎ　（Ｎ は、４で割って１余る数）　において、ｘ が偶数のときと同様に、

　ｘ²＋１＝２・Ａ・Ｂ　（Ａ は、４ｎ＋３型の素数で、Ｂ は整数）と分解されるものとする。

　このとき、Ａ・Ｂは、４で割って１余る数なので、Ｂは、４で割って３余る数でなければならな

い。よって、適当に、Ｂの素因数のうち ４ｎ＋３型の素数とＡとを入れ替えることにより、

Ａ≦Ｂ　としても一般性を失わない。

　ここで、　ｘ＜Ａ　とすると、　ｘ＋１≦Ａ≦Ｂ　で、

　　　ｘ²＋１＝２・Ａ・Ｂ≧２（ｘ＋１）²＝２ｘ²＋４ｘ＋２

となり、矛盾するので、　Ａ≦ｘ　である。

　このとき、　（ｘ−Ａ）²＋１＝ｘ²＋１−２ｘＡ＋Ａ²＝２・Ａ・Ｂ−２ｘＡ＋Ａ²＝Ａ（Ａ＋２Ｂ−２ｘ）

　ｘ−Ａ＜ｘ　で、帰納法の仮定より、（ｘ−Ａ）²＋１ の奇数の素因数は、４ｎ＋１型の素数の

みのはず。ところが、奇数 Ａ は、４ｎ＋３型の素数なので矛盾。

　したがって、３以上の奇数 ｘ に対して、命題は成立する。

　以上から、全ての自然数 ｘ に対して、ｘ²＋１ の奇数の素因数は、４ｎ＋１型の素数のみ

であることが示された。（証終）


　このことから、次が成り立つことは明らかだろう。

系　（偶数）²＋１　の素因数はすべて ４ｎ＋１型の素数である


　さらに、フェルマーにより見出され、オイラーにより証明された次の事実は有名であろう。

　　４ｎ＋１型素数は、２つの平方数の和として表せる

　平成２２年７月２５日付けで、当ＨＰの掲示板「出会いの泉」にＨＮ「ｋ．ｎｉｋａ」さんが、この
４ｎ＋１型素数についての話題を書き込まれた。

　次のようなものを考えてみました。

　すべての正の整数は高々３つの三角数で表す事が出来るが、特に

４ｎ＋１型素数の ｎ は　２つの三角数の和で表す事が出来る。

       ｎ＝ａ（ａ＋１）/ ２ ＋ ｂ（ｂ＋１）/ ２          (ただし、ａ、ｂ は非負の整数）

（注意）　便宜的に、０（０＋１）/２＝０　も三角数と呼ぶことにする。

（証明）　４n+１型素数は、２個の平方数の和で表すことができるので、Ｐ＝Ａ²＋Ｂ²　と置ける。

　　　両辺を２倍して右辺を式変形し、

　　　　　２Ｐ＝２Ａ²＋２Ｂ²＝（Ａ＋Ｂ)²＋(Ａ−Ｂ )²   （ただし、Ａ＞Ｂ）

　　Ｐ＝Ａ²＋Ｂ²　より、Ａ と Ｂ について、一方は偶数、他方は奇数になるので、

　　　　　Ａ＋Ｂ＝２ａ＋１　、Ａ−Ｂ＝２ｂ＋１　　（ただし、ａ、ｂ は非負の整数）

　とおける。このとき、

　　　　２Ｐ＝（２a＋１）²＋(２b＋１）²＝４a²＋４a＋１＋４ｂ²＋４ｂ＋１

　　よって、　Ｐ＝２ａ²＋２ａ＋２ｂ²＋２ｂ＋１

　　　　Ｐ＝４ｎ＋１　なので、　４ｎ＋１＝２ａ²＋２ａ＋２ｂ²＋２ｂ＋１

　　したがって、　ｎ＝ａ（ａ＋１）/２ ＋ ｂ（ｂ＋１）/２　　となる。　　（証終）

（コメント）　Ｐ＝Ａ²＋Ｂ²　から　Ｐ＝２ａ²＋２ａ＋２ｂ²＋２ｂ＋１　への精密化から得られる結
　　　　　　果ですね！ｋ．ｎｉｋａ さんに感謝します。

　さらに、ｋ．ｎｉｋａ さんは、次の問題も考えられた。

　４ｎ＋１型素数Ｐは２つの平方数の和で表す事が出来るが、その２つの平方数は、

　　　　Ｐ＝（ａ＋ｂ＋１）²＋（ａ−ｂ）²　　（ただし、ａ、ｂ は非負の整数で　ａ＞ｂ）

と表す事が出来る。

例　　５＝２²＋１²＝（１＋０＋１）²＋（１−０）²

　　　１３＝３²＋２²＝（２＋０＋１）²＋（２−０）²

　　　１７＝４²＋１²＝（２＋１＋１）²＋（２−１）²

　　　２９＝５²＋２²＝（３＋１＋１）²＋（３−１）2

（証明）　４ｎ＋１型素数 Ｐ＝Ａ²＋Ｂ²　において、

　　　　　２Ｐ＝２Ａ²＋２Ｂ²＝（Ａ＋Ｂ)²＋(Ａ−Ｂ )²   （ただし、Ａ＞Ｂ）

　　このとき、Ａ＋Ｂ＝２ａ＋１　、Ａ−Ｂ＝２ｂ＋１　　（ただし、ａ、ｂ は非負の整数）　とおける

　ので、２Ｐ＝（２a＋１）²＋(２b＋１）²＝｛（a＋ｂ＋１）＋（ａ−ｂ）｝²＋｛（a＋ｂ＋１）−（ａ−ｂ）｝²

ここで、　ａ＋ｂ＋１＝α　、ａ−ｂ＝β　と置くと、２Ｐ＝（α+β）²＋(α−β）²＝２α²＋２β²

となり、　Ｐ＝α²＋β²＝（ａ＋ｂ＋１）²＋（ａ−ｂ）²　と表すことができる。　　（証終）


　この問題について、らすかるさんが別証を与えられた。（平成２２年７月２５日付け）

（別証）　Ｐ＝Ａ²＋Ｂ²　（ただし、Ａ＞Ｂ）　において、

　　　　　　　ａ＝（Ａ＋Ｂ−１）/２　、　ｂ＝（Ａ−Ｂ−１）/２

　　　　とおくと、 ａ、ｂ は非負の整数で、

　　　　　　　Ａ＝ａ＋ｂ＋１　、　Ｂ＝ａ−ｂ

　　　　なので、　Ｐ＝（ａ＋ｂ＋１）²＋（ａ−ｂ）²　と表すことができる。　　（別証終）

　さらに、ｋ．ｎｉｋａ さんは、次の問題も考えられた。

　４ｎ＋１型素数　Ｐ＝Ａ²＋Ｂ²　（ただし、Ａ＞Ｂ）　において、斜辺の平方がＰで、残り

の２辺がＡ、Ｂである直角三角形を考えるとき、その直角三角形の面積Ｓは、

　　　　Ｓ＝ａ（ａ＋１）/２ − ｂ（ｂ＋１）/２　　（ただし、ａ、ｂ は非負の整数　ａ＞ｂ）

と、高々２つの三角数の差に書ける。

（証明）　Ｓ＝ＡＢ/２＝（ａ＋ｂ＋１）（ａ−ｂ）/２＝（ａ²＋ａ−ｂ²−ｂ）/２　より、

　　　　Ｓ＝ａ（ａ＋１）/２ − ｂ（ｂ＋１）/２　　（証終）

　これについて、らすかるさんからのアドバイスです。（平成２２年７月２５日付け）

　三角数　（０）、１、３、６、・・・　の定義から、その階差数列は自然数列　１、２、３、・・・　と
なる。よって、面積が整数であれば、必ず２つの三角数の差で表すことができる。

　また、ｋ．ｎｉｋａ さんは、４ｎ＋１型素数Ｐにおいて、２つの平方数を求める表を作られた。
　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　（平成２２年７月２５日付け）

　この表の使い方について、４ｎ＋１型素数　Ｐ＝Ｙ²＋Ｚ²＝１０９　を例として説明します。

（１）　まず、最初に ｎ を求める。　・・・・・　ｎ＝（１０９−１）/４＝２７

（２）　表のやや中央、赤い三角数のラインから下で、「２７」を探す。

（３）　「２７」を見つけたら、そこから垂線を下ろすと、赤いＹ列の交点の数字とぶつかる。

　　　この場合、「１０」が Ｙ の値となる。

（４）　次に、（２）で探した「２７」から、左斜め上４５度のラインで同じ「２７」を探す。

（５）　探したら、その数字の真上のＺ列との交点の数字「３」が Ｚ の値となる。

（確認）　１０２＋３２＝１００＋９＝１０９

　　　　ちなみに、この表の数字の意味はこちら。

（コメント）　何となく意味深な表ですね！ｋ．ｎｉｋａ さんに感謝します。


　ｋ．ｎｉｋａ さんが、「４ｎ＋１型奇数の ｎ」について新たな問題を提起された。
　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　（平成２２年８月１日付け）

　α²＋β²　という形の ４ｎ＋１型の奇数 Q を考える。 α＝２ｃ＋１、β＝２ｄ　 (ｃ、ｄ： 非負
整数）とする。そして、次のような操作を行う。

　（１）　ｎ を超えない最大の平方数を求め ｒ とする。（ｎ≧r)
　（２）　その数を ｎ から引き、その答えを ｅ とする。
　（３）　次に、ｅ を４倍して１を足し、その答えを ｆ とする。
　（４）　√ｆ を求める。

　この操作を、α の値を一定にし、β の ｄ値 を１ずつ増やす形で計算し、横に並べてみる。

α＝１の場合

　この場合、　√ｆ＝１　の値は　αの値　”１”　と一致している。

α＝３ の場合

　この場合も、√ｆ＝３　の値は、αの値　”３”　と一致している。

α＝５ の場合

　この場合も、√ｆ＝５　の値は、αの値　”５”　と一致している。

　もし、このまま上手くいけば、４ｎ＋１型奇数Ｑの中には、４ｎ＋１型素数Ｐも当然含まれる
ので、ｎ の値から最大平方を引くだけで、Y値が分かり、おのずと Z値も分かるのではない
かと期待してしまう。

α＝７ の場合

　この場合、前半２つが破綻してしまったが、後の√ｆ＝７の値は、αの値”７ ”と一致してい
る。因みに、破綻した ｎ＝３７ に対しては、最大平方より１つ下の平方数 ５² を引いてみる
と、３７−５²＝１２　から、　ｆ＝４９　となり、√ｆ＝７　で α値と一致する。

　ｎ＝２８に対しても同様に計算して、２８−4² ＝１２　で、√ｆ＝７　となる。

ここで、この破綻の仕組みについてご教授をお願いしたい。

　（１）　αの値と破綻する番号との因果関係
　（２）　一度破綻から回復すると、もう二度と破綻は起こらないのか
　（３）　４ｎ＋１型素数 Ｐ＝Y²＋Z²　のＰのみが与えられた場合、ｎ から、最大平方から何
　　　　番目の平方数を引けば、YとZが求まるか分かるアルゴリズム
　（４）　４ｎ＋１型素数Ｐが与えられた場合、ｎ から最大平方を引いてY、Zが求まる確率は？


　上記について、らすかるさんからコメントをいただきました。（平成２２年８月２日付け）

　（１）　αの値と破綻する番号との因果関係

　　　破綻するのは α²≧４β＋５ つまり β≦(α²−５)/４ の場合。

　　　α＜βに限定せず、β＝２{(番号)−１} とするのであれば、(番号)≦(α²＋３)/８

　　のとき破綻する。

　　　α＜βに限定する場合は、 (番号)=(β−α＋１)/２ なので、(番号)≦{(α−２)²−５}/８

　　のときに破綻する。

　　　なぜα＜βに限定している？

　（２）　一度破綻から回復すると、もう二度と破綻は起こらないのか

　　　上記により、起こらない。

　（３）　４ｎ＋１型素数 Ｐ＝Y²＋Z²　のＰのみが与えられた場合、ｎ から、最大平方から何
　　　　番目の平方数を引けば、YとZが求まるか分かるアルゴリズム

　　　それがわかったらすごいですね。

　（４）　４ｎ＋１型素数Ｐが与えられた場合、ｎ から最大平方を引いてY、Zが求まる確率は？

　　　α≪βである可能性が高いとは考えられないので、一般的には、α²≧４β＋５　を満
　　たし、確率は0になると思います。


　ＦＮさんも、４ｎ＋1型素数Ｐを α²＋β² の形で表すことについて、ExcelVBAで調べられま
した。（平成２２年８月２日付け）

　ｋ．ｎｉｋａ さんは、α を奇数、β を偶数としていますが、そうはしないで、α＞βとします。
100万以下で調べたら瞬時だったので1億以下で調べました。さすがにかなり時間がかかり
ましたが10分弱でした。むしろ書きだしたテキストファイルを読み込むのに苦労しました。適
当なエディターを持っていれば問題ないでしょうが、残念ながら持っていません。Excel2007
は100万行強まで読めますが、それでは足りません。ワードパッドで読みましたが非常に時
間がかかりました。

　１億の直近の10個の ４ｎ＋１ 型素数についての結果は次の通りでした。Ｒは、α² がＰ以
下の平方数のなかで大きい方から何番目かを表します。

  素数P　　　α　　β　　　R
99999517 , 7819 , 6234 , 2181
99999541 , 9890 , 1479 ,  110
99999589 , 9425 , 3342 ,  575
99999617 , 8129 , 5824 , 1871
99999677 , 7429 , 6694 , 2571
99999721 , 7795 , 6264 , 2205
99999773 , 8998 , 4363 , 1002
99999821 , 8765 , 4814 , 1235
99999941 , 7346 , 6785 , 2654
99999989 , 8458 , 5335 , 1542

　次に、１億以下を1000万で区切って、４ｎ＋1型素数のＲについて調べました。

      　区間　  　　　　素数の個数　R=1  R=1の確率　 Rの平均値
       　　1〜 9999999　 332180 14427 0.043431272　205.6074658
10000000〜19999999　302990　8815　0.029093369　384.2761543
20000000〜29999999　293609　7491　0.025513523　497.7437851
30000000〜39999999　287908　6704　0.023285216　589.6343624
40000000〜49999999　283765　6309　0.022233186　669.4256903
50000000〜59999999　280218　5872　0.020955114　739.0882242
60000000〜69999999　278048　5655　0.020338215　803.9461064
70000000〜79999999　275655　5295　0.019208794　863.0004897
80000000〜89999999　274020　5075　0.018520546　921.0241588
90000000〜99999999　272111　4934　0.018132306　972.00692

　　※素数の個数は、４ｎ＋1型の素数の個数、Ｒ＝１は、Ｒ＝１の個数。

　ｋ．ｎｉｋａ さんが問題にされていた確率は、Ｒ＝１の確率ですが、1億近くで２％弱です。５０
個に1個近くがＲ＝１だということで、かなり大きいとも言えます。因みに、1億以下で、Ｒ＝１
となるもっとも大きい素数は「99998497」です。「99999989」から小さい方へ数えて３５番目
です。1億と言えば普通にはかなり大きい数ですが、たかが10⁸だとも言えます。このぐらい
の桁数なら今のコンピュータならごく単純なプログラムでも処理できそうです。10¹⁵ぐらいに
なるとどうでしょうか。「1000000000000037」は、Ｍａｘｉｍａで調べると素数で、４ｎ＋１型です。
これを、α²＋β² の形で表すα、βを求めることはできるでしょうか。

　ｋ．ｎｉｋａ　さんより（平成２２年８月３日付け）

　どうもありがとうございます。αの値と破綻する番号との因果関係が不規則で、悩んでおり
ました。早速らすかるさんの計算式をもとに、さらに発展させていきたいと思います。
　究極的には、４ｎ＋１型素数Ｐの値のみから、Y²＋Z² を求める方法を考えたいと思います。
また、FNさん、どうもありがとうございます。大変時間のかかる作業、感謝いたします。このデ
ータはとても助かります。R＝１が５０個に１個とはとても高い確率だと思います。R＝１〜　が
それぞれ単純に同じ確率だとして、たとえば、ある ４ｎ＋１型素数Ｐが与えられたとして、R＝
１０まで計算していいとすると、５個に１個は、Y²＋Z² が求められるという事ですから、これは
すばらしいデータであると思います。あとは、「1000000000000037」を、Y²＋Z² に出来るよう
がんばります。

　らすかるさんが、より大きな素数について調べられた。（平成２２年８月３日付け）

10¹⁵〜10¹⁵+1000　の中で、

1000000000000037=26043586²+17936879²　、R=5579191
1000000000000241=31501204²+2770225²　　、R=121573
1000000000000249=30465032²+8478315²　　、R=1157745
1000000000000273=31508097²+2690692²　　、R=114680
1000000000000297=29592564²+11148101²　　、R=2030213
1000000000000357=27549121²+15525654²　　、R=4073656
1000000000000513=31506672²+2707327²　　、R=116105
1000000000000613=23604367²+21043618²　　、R=8018410
1000000000000741=28153375²+14400954²　　、R=3469402
1000000000000873=23018608²+21682797²　　、R=8604169
1000000000000921=29491515²+11412736²　　、R=2131262

10¹⁶〜10¹⁶+1000　の中で、

10000000000000061=86561206²+50071525²　　、R=13438795
10000000000000069=91785362²+39691905²　　、R=8214639
10000000000000453=99838442²+5682033²　　、R=161559
10000000000000481=95772809²+28767500²　　、R=4227192
10000000000000597=91441289²+40478274²　　、R=8558712
10000000000000613=96802783²+25084282²　　、R=3197218
10000000000000669=87688190²+48070587²　　、R=12311811
10000000000000753=98285017²+18440592²　　、R=1714984
10000000000000793=80486723²+59345492²　　、R=19513278
10000000000000861=99890394²+4680725²　　、R=109607
10000000000000897=71643416²+69765471²　　、R=28356585
10000000000000909=98554845²+16939378²　　、R=1445156
10000000000000949=91670018²+39957575²　　、R=8329983
10000000000000957=90930834²+41612299²　　、R=9069167

10¹⁷〜10¹⁷+1000　の中で、

100000000000000013=297277123²+107825378²　　、R=18950644
100000000000000021=270034886²+164563545²　　、R=46192881
100000000000000049=247146257²+197278300²　　、R=69081510
100000000000000081=291851400²+121748759²　　、R=24376367
100000000000000141=284313570²+138440579²　　、R=31914197
100000000000000181=316145770²+7200841²　　、R=81997
100000000000000337=311490196²+54533089²　　、R=4737571
100000000000000369=315693300²+18377713²　　、R=534467
100000000000000589=309204842²+66274925²　　、R=7022925
100000000000000609=315940728²+13470575²　　、R=287039
100000000000000669=284366725²+138331362²　　、R=31861042
100000000000000781=262872595²+175778266²　　、R=53355172
100000000000000817=284150096²+138775801²　　、R=32077671
100000000000000889=284524708²+138006125²　　、R=31703059

10¹⁸〜10¹⁸+1000　の中で、

1000000000000000009=1000000000²+3²　　、R=1
1000000000000000177=773571884²+633708561²　　、R=226428117
1000000000000000201=996077099²+88489620²　　、R=3922902
1000000000000000381=943385285²+331698966²　　、R=56614716
1000000000000000621=886461614²+462802125²　　、R=113538387
1000000000000000841=1000000000²+29²　　、R=1
1000000000000000861=907385115²+420300194²　　、R=92614886
1000000000000000877=983652806²+180075421²　　、R=16347195
1000000000000000913=977777983²+209643068²　　、R=22222018
1000000000000000997=997964566²+63770879²　　、R=2035435

　数が大きくなると、Ｒもどんどん大きくなりますが、10¹⁸〜10¹⁸+1000に、Ｒ＝１が２つもある
のは奇跡的な感じがします。

　理論値を概算してみました。

　Ｎが大きい値として、α²＋β² がＮ以下になる組合せは、α＞βとして約(π/８)N通りあり、
このうち、Ｒ＝１となる組合せは約(２/３)N^(3/4)通りあります。よってこの中に素数が均等
に出現すると仮定すると、Ｎ以下の素数で、Ｒ＝１となる確率は、

　　　{(２/３)N^(3/4)}/{(π/８)Ｎ}＝(１６)/(３πＮ^(1/4))≒２．４/Ｎ^(1/4)

となります。

　　N＜10000000 で 約0.043　　、　N＜100000000 で 約0.024

となりますので、FNさんが調査された値（約２％弱）とほぼ合っています。

　「Ｒ＝１が５０個に１個とはとても高い確率」について、確かにＮの大きさの割にはそのよう
に思えます。４乗根ですからね。しかし、N→∞　のとき、(確率)→０　ですから、例えば、５０
桁のような大きい数になると「ほぼ皆無」と言えますね。

　らすかるさんの結果に対して、ＦＮさんからのコメントです。（平成２２年８月３日付け）

　10¹⁵ ぐらいで大変かと思ったのですが、その辺は軽く超えて、10¹⁸ でもできますか。10²⁵
ぐらいになればさすがにきついでしょうか。Excel では15桁までなので、1000000000000037
と入れると1000000000000030になってしまいます。いつまでも ExcelVBA ばかりではだめで
すね。今日VisualC++をInstallしました。使う気になるまで何日もかかりそうですが。

　４ｎ＋１型の素数Ｐは、２つの自然数の平方の和に書けるという定理の証明は読みました。
オイラーによるものと本質的に同じだそうですが、その証明は実際に２つの平方の和に書く
手順を与えます。しかし相当に面倒です。だから普通はＰから平方数を引いていって、平方
数になるものを捜すのが実際的なようです。理論値というあたりの話は全くわかりません。

　らすかるさんが用いられたプログラムの概要です。（平成２２年８月４日付け）

　４ｎ＋１＝Ｎ　に対して、α＝[√N]、β＝[√(N-[√N]²)]　として、

　　α²＋β²　がNより小さければ、βに１を足す
　　α²＋β²　がNより大きければ、αから１を引く

以上２行を繰返し、α²＋β²＝Ｎ　となったら終わり。

　乗算をしなくて良いように若干工夫はしていますが、基本は上記だけです。

　10¹⁵　あたりは一瞬ですが、10¹⁸　あたりでは時間がかかるもの（Ｒが大きいもの）は２０秒
ほどかかりました。

　また、らすかるさんは、ＦＮさんの「１０²⁵ ぐらいになれば、さすがにきついでしょうか？」とい
う問いかけに答えて、プログラムを作り直され、次の計算結果を得られた。

　新しいアルゴリズムでは、最後の１０²⁰⁰+357 でも１秒もかかりません。

　　　参考ページ: http://ror.hj.to/r/blog/category/1/16/1

10²⁵〜10²⁵+1000　の中で、

10²⁵+13=2997606971002²+1007150657747²
10²⁵+349=2966451518443²+1095520601690²
10²⁵+513=3119980365017²+515482804668²
10²⁵+561=2339298527015²+2127835144344²
10²⁵+609=2703063575820²+1641172539703²
10²⁵+657=3093178628856²+657454158089²
10²⁵+937=3152592111584²+247311499891²

10⁵⁰〜10⁵⁰+1000　の中で、

10⁵⁰+577=7611065343808354245450401²+6486268906873921642245424²
10⁵⁰+709=9250603407870174268436578²+3798201757450585752871775²
10⁵⁰+889=8360474083778077443654955²+5486572071382560131393192²
10⁵⁰+897=9182579133974339811099231²+3959828335710913551362744²
10⁵⁰+961=10000000000000000000000000²+31²

10¹⁰⁰〜10¹⁰⁰+1000　の中で、

10¹⁰⁰+949
=99697921470138519447541656418848509184628524016382²
　　　　　　　　　+7766881905507050845172598218029833369440123277895²

10²⁰⁰〜10²⁰⁰+1000　の中で、

10²⁰⁰+357
=87342995146970964150923244222621726434876511735966
　10058044551075223059534652741558693062858574717874²
　　　　　　　　+48694981248134869822313779563554735739904013676676
　　　　　　　　　05534476012902371537410710789858816644285044662191²

　ＦＮさんからのコメントです。（平成２２年８月４日付け）

　１０²⁰⁰ あたりでも瞬時にできてしまうとは全く驚くばかりです。参考ページを読ませていた
だきました。

　中心は「4ｎ＋１型の素数が２つの平方数の和である」の証明(本質的にオイラーによる)で
した。やっぱり、オイラーは偉い。オイラーの証明はベストに近い見つけ方でもあったようで
す。A²＋B²＝ＭＰ　から、ａ²＋ｂ²＝ｍｐ　を作るとき、ｍ≦Ｍ/２　ですが、例としてあがって
いるケースでは、おおかた１桁落ちてますね。そこまでいかなくても２回で１桁ぐらい落ちそ
う。この部分は時間に関しては問題なさそうです。

次に、ｘ²≡−１　を解くことですが、確かに原子根がわかっていれば問題ないです。(と言っ
ても仮に２が原子根だとして２^(ｐ-1)/2を計算するプログラムを私が書けるかと言えば書けま
せん。ｐが200桁ぐらいあるとして)では、原子根をどうやって求めるかですが、ここに書いて
あるのは、ｐ−１を素因数分解してます。200桁ぐらいある数の素因数分解は一般には無茶
苦茶な時間がかかると聞いてます。試しにMaximaでfactor(10^200+356);を実行したら黙っ
てしまいました。しかし原子根はいっぱいあるから、２、３、５、・・・を (ｐ-1)/2乗していって、
−１になるものを捜せばかなり早くみつかるのでしょう。

ということで理屈としてはできそうです。そのことと実際にプログラムを書くこととはかなり距
離がありますが。

　上の参考ページに参考としてあげてある次のページに、「４ｎ＋１型の素数が２つの平方数
の和である」の簡明な証明が書いてあります。証明だけであればこんなに簡単にできるのに
驚きました。ただ実際に２つの平方数をみつけるのには役に立たないです。やはりオイラーに
よる証明でないと。

　さらに一文証明(one-sentence proof)なるものが載ってますが、実際には確認すべきことが
非常に多く one-page proof ぐらいになりそうです。でもおもしろそうです。なお、明らかですが
一文証明の６行目のif x<2yはif x>2yの間違いです。

らすかるさんからのコメントです。（平成２２年８月４日付け）

　私は最初そのページを見て、「桁数が多くなると素因数分解が困難になるからダメだ」と思
ったのですが、実際に試しているうちに、素因数分解は必要なく、また原始根を見つける必
要もないことに気付きました。

　FNさんのおっしゃる通りで、２、３、５、・・・を　(p-1)/2乗していって−１になるものを探せば
良いだけです。(p-1)/2乗が−１になるものは原始根とは限りませんが、Ａ＝ａ^(p-1)/4 とすれ
ば、Ａ²＋１²＝ｋｐ　となるのですから原始根でなくても何も問題ありません。実際に手作業
で百個以上の素数で調べた中では、(p-1)/2乗して−１になるものが17までで必ず見つかっ
ていました。

　ｋ．ｎｉｋａ さんからのコメントです。（平成２２年８月５日付け）

　らすかるさんからの質問：「なぜ、α＜β に限定しているのですか？」について、

（α＞β）において、β＝４の場合、ｓ ＝（ ｎ を（√ｆ ＝α）に変えるための平方数 ）とする。

　α＞βの場合には、Nから最大平方を引くのではなくて、√ｆ＝αにするための ｎ から引く、
あるひとつの平方数を求めなくてはいけません。よって、この場合は、R＝１の確率は非常に
低くなってしまいます。もし、（α＞β）の集合と（α＜β）の集合のそれぞれのR＝１を求める
とぜんぜん違う値になると思います。

　ＦＮさんからのコメントです。（平成２２年８月５日付け）

　一度原子根を求めたいと思っていたので、(p-1)/2乗して-1になるものとを調べてみました。
100万以下の素数Pについて2、3、5、･･･と調べて最初に見つかる原子根Gと、最初に見つか
る(p-1)/2乗が-1になるものFを調べてみました。ただし、４ｎ＋１型ではなくすべての奇素数
について調べました。

　F＝２　となるものが半分ぐらいあります。最大のFは４３でした。Fについては素数だけ調べ
れば十分です。小さい方から何番目の素数がFになるかの平均値は、１．９９１ほどでした。
４ん＋１型だけで調べた方がいいですがほとんど差はないでしょう。平均して２つ目でみつか
るわけだから全く問題なしです。参考までに、１９までについて、G＝A、F＝A の個数を書い
ておきます。

　また、「なぜα＜βに限定しているのですか？」について、「α＜βに限定してはまずいの
ではないですか」という意味だと思います。他の限定がなければα＜βとしていいですが、
ｋ．ｎｉｋａ さんは、α、β の偶奇を決めておられますから、その上で、α＜β とすると半分ぐ
らいが考慮の対象外となってしまいます。あとで、α＞β のケースは考えようということでし
ょうか。

　ｋ．ｎｉｋａ さんからのコメントです。（平成２２年８月６付け）

　あとで、α＞βのケースは考えようと考えております。α＞β のケースと α＜β のケース
では、α の求め方がまったく違うので別々に考えた方がいいのでは？と思うのです。
α＞βの集合について１億ぐらいまではR＝１〜１０までの確率は二桁（％）いくのではないか
と思っております。

　らすかるさんからのコメントです。（平成２２年８月６付け）

　G＝A、F＝Aの個数は私も気になっていて、自分でも調べようと思っていました。今までの
方法でも十分効率が良いのですが、効率をさらに上げる方法に気づきました。

　もし、(p-1)/4 が偶数の場合、ａ^(p-1)/4 が -1 でも良いわけですし、(p-1)/8 も偶数なら、
ａ^(p-1)/8 が -1 でも良いわけです。

　ａ^N　(ｍｏｄ　ｐ)を計算する場合、もともと　(p-1)/2^k＝q　（奇数）で、ａ^ｑ　(ｍｏｄ　ｐ)を計算し
てから2乗を繰り返していますので、上記の値は計算途中で既に得ています。

　よって、次のように処理すると効率良く計算できます。

　・ａ^ｑ を計算し、値が 1 か -1 なら、ａ を変更して繰り返します。

　・ａ^ｑ が 1 でも -1 でもなければ、これを順次２乗していけば、いつか -1 となり、その時の
　　２乗する前の値がAとなります。

# もし、-1 とならずに 1 になったり、(p-1)/2乗まで計算しても -1 にならなかった場合は、
　元の数は合成数です。

　この改善策の効果は以下のようになりました。
（対象は10億以下の4k+1型の数で、上記の判断で合成数と判断されなかった合成数（＝偶
然-1になったもの）もわずかですが含みます。）

　かなり改善されていることがわかると思います。改善前は平均２個目でしたが、改善後は
平均１．２個目です。改善後のプログラムで10⁵⁰⁰+1189を平方和に分けたら、実行時間は1
秒でした（改善前は１．６秒）。因みに、10⁵⁰⁰+1189は、10⁵⁰⁰より大きい２番目の素数です。
１番目は10^500+961なのですが、これは分解結果が、(10²⁵⁰)²+31² となってあまりにも面
白くないので次のを探しました。ここら辺になると、素数を探すのも大変ですね。


　ｋ．ｎｉｋａ さんからの新しい問題提起です。（平成２２年８月１７日付け）

　またひとつ　ご教授頂きたい事がございます。違う方法を考えてみました。

　２つの４ｎ＋１型素数を、Ｐ_n、Ｐ_m とし、　Ｐ_n＝Ａ²＋Ｂ²　、Ｐ_m＝Ｃ²＋Ｄ²　とする。

　Ｐ_n・Ｐ_m＝Q　とし、Q の値のみ与えられているものとする

　Qを平方和に分けるアルゴリズムを用いて、２種類の平方和　Q₁、Q₂ に分ける。

（ブラーマグプタの二平方恒等式- Wikipedia）

　　Q₁＝（AC - BD）²＋（AD＋BC）²　、Q₂＝（AC＋BD）²＋（AD - BC）²

このとき

　（１）　この２種類以外平方和は存在しないのか。

　（２）　この２種類の平方和は、ブラーマグプタの二平方恒等式の内容と一致していると見
　　　て良いのか。

　（３）　探し出した４つの平方数から、どれが（AC - BD）で（AD＋BC）であるとか調べる方
　　　法はあるのか。（これは分かればラッキーという感じです。）

　以上３点になります、これが肯定的に分かれば、うまくすれば素因数分解の桁数を減らせ
そうな気がするのですが・・・どうぞ宜しくお願いいたします。


　この問いかけに対して、らすかるさんからのコメントです。（平成２２年８月１８日付け）

　（１）　この２種類以外平方和は存在しないのか。

　　　４ｎ＋１型の２素数の積であれば、実質的に２通りしか存在しないはずです。

　（２）　この２種類の平方和は、ブラーマグプタの二平方恒等式の内容と一致していると見
　　　て良いのか。

　　　２通りですから、一致していると思います。

　（３）　探し出した４つの平方数から、どれが（AC - BD）で（AD＋BC）であるとか調べる方
　　　法はあるのか。

　　　(AC-BD)と(AC+BD)の偶奇は同じ、(AD-BC)と(AD+BC)の偶奇は同じで、

　　　　　(AC-BD)＜(AC+BD)　、 (AD-BC)＜(AD+BC)

　　 ですから、CとDの入れ替えを除けばただちにわかりますね。

　ただし、問題は大きな合成数の場合は平方和に分けるのが困難であるという点です。

　例の平方和に分けるアルゴリズムは、素数の場合は、−１の平方根が容易に見つかりま
すので巨大な数でも高速に処理できますが、合成数の場合はそうはいきません。もし、−１
の平方根が容易に見つかるのであれば、平方和に分けなくてもただちに素因数分解できま
す。

　このことについて、ｋ．ｎｉｋａ さんからのコメントです。（平成２２年８月１９日付け）

　平方和に分けるには、−１の平方根を求める必要があり、−１の平方根が見つかれば、
平方和に分ける必要はないとは困りました

　では、そういう難しい議論はとりあえず抜きにして、２つの４ｎ＋１型素数の積を、方法は問
わず２種類の平方和に分ける事が出来たという流れで、私の気がついた事について書かせ
て頂きたいと思います。

　まず、２つの４ｎ＋１型素数を、Ｐ_n、Ｐ_m とし、　Ｐ_n＝Ａ²＋Ｂ²　、Ｐ_m＝Ｃ²＋Ｄ²　とする。

　次に、　Ｐ_n・Ｐ_m＝Q　とし、Q の値のみ与えられたとする

　Qを２種類の平方和に分ける。

　　Q₁＝（AC - BD）²＋（AD＋BC）²　、Q₂＝（AC＋BD）²＋（AD - BC）²

　次に、以下の計算をする。

　　　（（AC - BD)＋(AC＋BD））/2＝AC

　　　（（AD - BC)＋(AD＋BC））/2＝AD

　　　（（AD＋BC)-(AD - BC））/2＝BC

　　　（（AC＋BD)-(AC - BD））/2＝BD

　　　AC＋AD＋BC＋BD＝（ A＋B )( C＋D ）

　ここまでの目的は、　（Ａ²＋Ｂ²）（Ｃ²＋Ｄ²）を（Ａ＋Ｂ )( Ｃ＋Ｄ ）　の形に変える事です。

　次に、（Ａ＋Ｂ )( Ｃ＋Ｄ ）を既存の素因数分解プログラムで解く。

　　　（Ａ＋Ｂ )( Ｃ＋Ｄ ）＝ F^a×G^ｂ×H^ｃ×I^ｄ×・・・

　次に、その素因子をいろいろ組み合わせて、以下の式に入れて割り切るまで試す。

　たとえば、　（ ４＋５ )( ７＋８ ）＝１３５＝３×３×３×５＝（３×３)(３×５）とする。

　まず、（３×５）を（ Ｃ＋Ｄ )であると仮定して、

　　　｛（ＡＣ ＋ ＡD）²＋（ＢＣ＋ＢＤ）²｝/（Ｃ＋Ｄ）²＝Ａ²＋Ｂ²

　これで、Qを割ってみて割り切れれば正解になります。（以上）

　この計算で、うまくいけば、たとえば２００桁の計算を１００桁強ぐらいにまで落とせるので
はと思うのですが．．．。なにか計算ミス、見落とし等あればご教授頂けますよう、宜しくお
願いいたします。

　らすかるさんからのコメントです。（平成２２年８月１９日付け）

　間違いはないと思いますが、ＡＣ、ＡＤ、ＢＣ、ＢＤが出れば、

　　Ａは、ＡＣとADの最大公約数　、Ｂは、ＢＣとＢＤの最大公約数、

　　Ｃは、ＡＣとＢＣの最大公約数　、Ｄは、ＡＤとＢＤの最大公約数

ですから、さらに素因数分解する必要はないですね。

実例　　790037は、790037=71²+886²　、790037=626²+631²　の２通りに書ける。

　偶奇が同じ組の和と差の半分を求めて、

　(631+71)/2=351　、 (631-71)/2=280　、 (886+626)/2=756　、 (886-626)/2=130

　GCD(351,756)=27　、GCD(280,130)=10　、GCD(351,130)=13　、GCD(280,756)=28

　　∴　790037=(27²+10²)*(13²+28²)=829*953　のように求められる。

　因みに、　Ｘ²≡−１ (ｍｏｄ 790037) となる790037/2以下の２数は、157756 と 255915です
が、これらがわかっていれば、

　　GCD(790037,255915+157756)=829　、GCD(790037,255915-157756)=953

のように素因数が求められます。


（追記）　平成２２年９月２０日付け

　この話題については、ｋ．ｎｉｋａさんのＨＰサイト「４ｎ＋１型素数の部屋」をご覧ください。
分かり易く説明されています。



　　以下、工事中