Prime Numbers in a sequence (1+2*{(1+sqrt(-7))/2}^k)*(1+2*{(1-sqrt(-7))/2}^k)

Prime Numbers in a sequence (1+2{(1+sqrt(-7))/2}^k)(1+2*{(1-sqrt(-7))/2}^k)

[2012.07.14]J_k=(1+2{(1+sqrt(-7))/2}^k)(1+2*{(1-sqrt(-7))/2}^k)が素数になるkを探す

J_kが素数になるkを探すプログラム(pari/gp)

■α=(1+sqrt(-7))/2とする。正整数kに対して、j_k, J_kを以下のように定義する。
     j_k = 1+2*α^k
     J_k = j_k*conj(j_k) = 1 + 2*(α^k+conj(α)^k) + 2^k+2

ただし、conj(z)はzの共役複素数である。
また、α, j_kは2次の代数的整数である。

■J_kは、以下の線型漸化式を満たすことが分かる。
     J_k+4=4*J_k+3-7*J_k+2+8*J_k+1-4*J_k
     J₁=J₂=1, J₃=23, J₄=67
よって、各k>0に対して、J_kは有理整数である。

gp > read("jk.gp")
gp > JJJ(1)
time = 0 ms.
%1 = 11
gp > JJJ(2)
time = 0 ms.
%2 = 11
gp > JJJ(3)
time = 0 ms.
%3 = 23
gp > JJJ(4)
time = 0 ms.
%4 = 67

■参考文献[1]の方法に従って、数列{J_k}中の素数を探す。
[1]では、k<=700,000の範囲まで調べて、50個の素数を確定している。

■{J_k mod 3}は周期8, {J_k mod 5}は周期24である。

gp > for(k=1,16,print1(JJ(k)%3);print1(", "));print("...")
2, 2, 2, 1, 1, 2, 1, 0, 2, 2, 2, 1, 1, 2, 1, 0, ...
time = 0 ms.
gp > for(k=1,48,print1(JJ(k)%5);print1(", "));print("...")
1, 1, 3, 2, 1, 0, 2, 3, 4, 1, 2, 1, 2, 3, 3, 3, 2, 4, 4, 2, 4, 3, 4, 4, 1, 1, 3, 2, 1, 0, 2, 3, 4, 1, 2, 1, 2, 3, 3, 3, 2, 4, 4, 2, 4, 3, 4, 4, ...
time = 2 ms.

さらに、k≡0 (mod 8)のとき、J_kは3で割り切れ、k≡6 (mod 24)のとき、J_kは5で割り切れる。
よって、k≡0 (mod 8)またはk≡6 (mod 24)のとき、J_kが合成数であることが分かる。

■虚数乗法(CM)を持つ楕円曲線E
E: y² = x³ - 35*x + 98
の(-a)-twist
E_a: y² = x³ - 35*a²*x - 98*a³
で、E_a(Q)が正のrankを持つものを考える。ただし、aは平方因子を持たない整数とする。
楕円曲線E_a上の有理点P_aを以下のように定義する。

k	a	P_a
k≡0,2 (mod 3)	-1	(1,8)
k≡4,7,13,22 (mod 24)	-5	(15,50)
k≡10 (mod 24)	-6	(21,63)
k≡1,19 (mod 24)	-17	(81,440)
k≡25,43,49,67 (mod 72)	-111	(-633,12384)

gp > ec(a)=ellinit([0,0,0,-35*a^2,-98*a^3])
time = 0 ms.
%21 = (a)->ellinit([0,0,0,-35*a^2,-98*a^3])
gp > e1=ec(-1)
time = -2 ms.
%22 = [0, 0, 0, -35, 98, 0, -70, 392, -1225, 1680, -84672, -1404928, -3375, [-7.0000000000000000000000000000000000000, 3.5000000000000000000000000000000000000 - 1.3228756555322952952508078768196302129*I, 3.5000000000000000000000000000000000000 + 1.3228756555322952952508078768196302129*I]~, 2.5575309899160994790492257969408742850, 1.2787654949580497395246128984704371425 - 0.48332792640420288668326920975745343276*I, 3.9423885925549265795997104583799551347 + 2.350988701644575016 E-38*I, 1.9711942962774632897998552291899775673 + 1.7116972661997117051282981581331545224*I, 1.2361261500706366841245492420460944111]
gp > ellisoncurve(e1,[1,8])
time = 0 ms.
%23 = 1
gp > e2=ec(-5)
time = 0 ms.
%24 = [0, 0, 0, -875, 12250, 0, -1750, 49000, -765625, 42000, -10584000, -21952000000, -3375, [-35.000000000000000000000000000000000000, 17.500000000000000000000000000000000000 - 6.6143782776614764762540393840981510643*I, 17.500000000000000000000000000000000000 + 6.6143782776614764762540393840981510643*I]~, 1.1437626296029455197711153946718935103, 0.57188131480147275988555769733594675516 - 0.21615081977276263014541000326271581939*I, 8.8154488866725371360290702569558999558 + 0.E-37*I, 4.4077244433362685680145351284779499779 + 3.8274714441231084874878206264715746989*I, 0.24722523001412733682490984840921888222]
gp > ellisoncurve(e2,[15,50])
time = 0 ms.
%25 = 1
gp > e3=ec(-6)
time = 0 ms.
%26 = [0, 0, 0, -1260, 21168, 0, -2520, 84672, -1587600, 60480, -18289152, -65548320768, -3375, [-42.000000000000000000000000000000000000, 21.000000000000000000000000000000000000 - 7.9372539331937717715048472609177812771*I, 21.000000000000000000000000000000000000 + 7.9372539331937717715048472609177812771*I]~, 1.0441076544415988952242884345177083639, 0.52205382722079944761214421725885418196 - 0.19731779968795962713738947184403518978*I, 9.6568404195287026285644720787256250348 + 4.701977403289150032 E-38*I, 4.8284202097643513142822360393628125174 + 4.1927848963062009547008756449741872589*I, 0.20602102501177278068742487367434906852]
gp > ellisoncurve(e3,[21,63])
time = 0 ms.
%27 = 1
gp > e4=ec(-17)
time = 0 ms.
%28 = [0, 0, 0, -10115, 481474, 0, -20230, 1925896, -102313225, 485520, -415993536, -33911546540032, -3375, [-119.00000000000000000000000000000000000, 59.500000000000000000000000000000000000 - 22.488886144049020019263733905933713619*I, 59.500000000000000000000000000000000000 + 22.488886144049020019263733905933713619*I]~, 0.62029237718909259044378265679449454555, 0.31014618859454629522189132839724727277 - 0.11722424072795809049952274973249173981*I, 16.254884584334108807496706549132357239 - 9.403954806578300064 E-38*I, 8.1274422921670544037483532745661786193 + 7.0575086276224015793615455628075462398*I, 0.072713302945331569654385249532123200654]
gp > ellisoncurve(e4,[17,440])
time = 0 ms.
%29 = 0
gp > ellisoncurve(e4,[81,440])
time = 0 ms.
%30 = 1
gp > e5=ec(-111)
time = 0 ms.
%31 = [0, 0, 0, -431235, 134027838, 0, -862470, 536111352, -185963625225, 20699280, -115800052032, -2627797775938449408, -3375, [-777.00000000000000000000000000000000000, 388.50000000000000000000000000000000000 - 146.83919776408477777283967432697895363*I, 388.50000000000000000000000000000000000 + 146.83919776408477777283967432697895363*I]~, 0.24275009884636697002470432710773281988, 0.12137504942318348501235216355386640994 - 0.045875456591702455054197972254420731715*I, 41.535641170355139392810338777228467780 - 1.880790961315660013 E-37*I, 20.767820585177569696405169388614233890 + 18.033849726384766163044246890886959819*I, 0.011136271622257988145266209387802652352]
gp > ellisoncurve(e5,[-633,12384])
time = 0 ms.
%32 = 1
gp > ellpow(e1,[1,8],2)
time = 0 ms.
%33 = [2, -6]
gp > ellpow(e1,[1,8],3)
time = 0 ms.
%34 = [193, 2680]
gp > ellpow(e1,[1,8],4)
time = 0 ms.
%35 = [-47/144, 18073/1728]
gp > ellpow(e2,[15,50],2)
time = 0 ms.
%36 = [-26, -132]
gp > ellpow(e2,[15,50],3)
time = 0 ms.
%37 = [51615/1681, -8250850/68921]
gp > ellpow(e3,[21,63],2)
time = 0 ms.
%38 = [-167/4, -253/8]
gp > ellpow(e4,[81,440],2)
time = 0 ms.
%39 = [-132446/3025, 152522658/166375]
gp > ellpow(e5,[-633,12384],2)
time = 0 ms.
%40 = [66107953/29584, -517135087657/5088448]

■k!≡0 (mod 8)かつk!≡6 (mod 24)のとき、kに依存して、点P_a∈E_a(Q)を上記のようにとると、以下は同値になる。
(1)2^k+1P_a≡O mod J_k,
かつ、2^kP_a mod J_kの射影座標[x_k:y_k:z_k]に対して、gcd(z_k,J_k)=1.
(2)J_kが素数である。

よって、J_kが素数かどうかを楕円曲線E_a mod J_k上の有理点P_aの位数が2^k+1かどうかで判定することができる。
しかも、J_kが素数のとき、J_k+1, J_k-1は(2以外の)自明な素因数を持たない。

■Pari/GPでプログラムを作成して、k<=136600の範囲で素数となるJ_kを調べると、以下のように見つかった。

gp > check2(2,136600)
J_{2} is prime
J_{3} is prime
J_{4} is prime
J_{5} is prime
J_{7} is prime
J_{9} is prime
J_{10} is prime
J_{17} is prime
J_{18} is prime
J_{28} is prime
J_{38} is prime
J_{49} is prime
J_{53} is prime
J_{60} is prime
J_{63} is prime
J_{65} is prime
J_{77} is prime
J_{84} is prime
J_{87} is prime
J_{100} is prime
J_{109} is prime
J_{147} is prime
J_{170} is prime
J_{213} is prime
J_{235} is prime
J_{287} is prime
J_{319} is prime
J_{375} is prime
J_{467} is prime
J_{489} is prime
J_{494} is prime
J_{543} is prime
J_{643} is prime
J_{684} is prime
J_{725} is prime
J_{1129} is prime
J_{1428} is prime
J_{2259} is prime
J_{2734} is prime
J_{2828} is prime
J_{3148} is prime
J_{3230} is prime
J_{3779} is prime
J_{5537} is prime
J_{5759} is prime
J_{7069} is prime
J_{7189} is prime
J_{7540} is prime
J_{7729} is prime
J_{9247} is prime
J_{10484} is prime
J_{15795} is prime
J_{17807} is prime
J_{18445} is prime
J_{19318} is prime
J_{26207} is prime
J_{27140} is prime
J_{31324} is prime
J_{36397} is prime
J_{47294} is prime
J_{53849} is prime
J_{83578} is prime
J_{114730} is prime
J_{132269} is prime
J_{136539} is prime

ただし、ここまで計算するには、かなり時間が掛かるので、kの範囲を分割して、少しづつ検索した。

■J₁₃₆₅₃₉は、このような(41103桁の)素数である。
この方法で、J₈₃₅₇₈, J₁₃₆₅₃₉が素数であることを確認するには、それぞれ約13分、約41分掛かった。

gp > is_prime_J(136539)
time = 41min, 39,848 ms.
%1 = 1
gp > is_prime_J(5759)
time = 835 ms.
%2 = 1
gp > is_prime_J(7069)
time = 1,382 ms.
%3 = 1
gp > is_prime_J(83578)
time = 13min, 23,746 ms.
%4 = 1
gp > log(JJ(83578))/log(10)
time = 948 ms.
%5 = 25160.087037595548295964296821073127177
gp > log(JJ(5759))/log(10)
time = 11 ms.
%6 = 1734.2338050201956656263497725078043272
gp > log(JJ(7069))/log(10)
time = 21 ms.
%7 = 2128.5830993400110313563477245968901923
gp > is_prime_J(53849)
time = 4min, 20,469 ms.
%8 = 1
gp > log(JJ(53849))/log(10)
time = 391 ms.
%9 = 16210.766296501051343455053219808673985

■MorainのECPP V3.4.1で、1735桁の整数J₅₇₅₉が素数であることを、このように確認できた。
しかし、ECPPでは、2129桁の整数J₇₀₆₉が素数であることを確認できなかった(quick compositeness test実行中に、segmentation faultを起こす)。

[2012.08.16追記]
44447桁の整数J₁₄₇₆₄₇が素数であることを確認できた。

[参考文献]

[1]Alexander Abatzoglou, Alice Silverberg, Andrew V. Sutherland, Angela Wong, "Deterministic Elliptic Curve Primerlity Proving for Special Sequence of Numbers", 2012.
[2]Carl Pomerance, "Very Short Primality Proofs", Math. Comp., Vol.38(177), p315-322, 1987.

Last Update: 2012.08.16

H.Nakao