1だけならぶ素数

このワークシートはMath by Codeの一部です。 <レプ・ユニット素数ってなんだっけ?> レプ・ユニット数は、反復する1、つまり、repeatedUnitsの略した言い方。 レプ・ユニット素数はどのくらいあるでしょう。 111=999/9=(1000-1)/9=(10^3-1)/9と表せるように、 一般に(10^n-1)/9で111...1のように1だけが並ぶ整数が作れます。 nが9973以下では、n=2,19,23,317,1031の5個だけだといわれています。 実験してみましょう。 juliaで素数判定関数をまた使います。次に、juliaでn=9999までの整数に対して レプ・ユニット数(10^n-1)/9に素数filterをかけよう。 [IN]julia #====================== N=[x for x in 2:30] Rp=map(n->(n,((10^n-1)÷9)),N) println(N) println(Rp) #====================== [OUT] [2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30] [(2, 11), (3, 111), (4, 1111), (5, 11111), (6, 111111), (7, 1111111), (8, 11111111), (9, 111111111), (10, 1111111111), (11, 11111111111), (12, 111111111111), (13, 1111111111111), (14, 11111111111111), (15, 111111111111111), (16, 1111111111111111), (17, 11111111111111111), (18, 111111111111111111), (19, -938527119301061290), (20, 862919959050249102), (21, 430646668853801415), (22, 207190227713669347), (23, 22264046724521073), (24, 222640467245210737), (25, 176766442039934975), (26, -281973810012822641), (27, -770099869716054016), (28, 497554224488149447), (29, 876265784057149667), (30, 564104918922807068)] <参考> (・素数リストの作り方、素数判定についてはこちら、  ・リストにfilterをかける方法についてはこちら) 残念なら桁溢れが19けたで起きている。 juliaの整数intは19ケタが限界なので、 BigInt()を使って、ケタの限界をこえるように改造しよう。 1が3の倍数個あると3の倍数になるからとばそう。 1が2より大きい偶数個(2n)あると、左右にわけた1がn個ならんだ数の倍数になるからとばそう。 いやいや、 それだけではない。 また、m,nに対するレプ・ユユニット数Rm,Rnがあるとき、mがnの倍数ならば、RmはRnの倍数になる。 すると、mが合成数ならば、Rmはmの因子nのRnを因数にもつから、合成数になる。 レプユニット素数Rpを探すには、pが素数であることは必要だが、十分でなない。 これって、メルセンヌ素数と同様に、mが素数であることは必要だが十分ではないときと同じだね。 質問:m,nに対するレプ・ユユニット数Rm,Rnがあるとき、mがnの倍数ならば、RmはRnの倍数になる。これはどうやって証明しますか。 レプユニット素数の候補数をならべよう。 Pythonの場合はJuliaより無名関数やmap,filterでリストを作ると、記述が長めになる。 だから、約数個数を求める関数divCountを作って、1行に詰め込みすぎないようにしよう。 2から100までの整数のリストMで、dicCountが2個になる整数Pが素数のリストとなる。 素数個の1を並べることで、レプ・ユニット素数の候補ができる。 #[IN]Python #================================= M=list(range(2,100)) divCount=lambda n:len(list(filter(lambda m:n % m==0, list(range(1,n+1))))) P=[x for x in M if divCount(x)==2] Rp=list(map(lambda n:(n,((10**n-1)//9)),P)) print(P) print(Rp) #================================== Juliaでは、divcountを関数にしなくても、filterの中の無名関数で簡単にかける。 #[IN]Julia #================================= P=filter(n->length(filter( m -> n % m==0,1:n))==2, 2:100)  #べき乗のBigInt()は底につける。 Rp=map(n->(n,(BigInt(10)^n-1)÷BigInt(9)), P) println(P) println(Rp) #================================= [OUT] [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97] [(2, 11), (3, 111), (5, 11111), (7, 1111111), (11, 11111111111), (13, 1111111111111), (17, 11111111111111111), (19, 1111111111111111111), (23, 11111111111111111111111), (29, 11111111111111111111111111111), (31, 1111111111111111111111111111111), (37, 1111111111111111111111111111111111111), (41, 11111111111111111111111111111111111111111), (43, 1111111111111111111111111111111111111111111), (47, 11111111111111111111111111111111111111111111111), (53, 11111111111111111111111111111111111111111111111111111), (59, 11111111111111111111111111111111111111111111111111111111111), (61, 1111111111111111111111111111111111111111111111111111111111111), (67, 1111111111111111111111111111111111111111111111111111111111111111111), (71, 11111111111111111111111111111111111111111111111111111111111111111111111), (73, 1111111111111111111111111111111111111111111111111111111111111111111111111), (79, 1111111111111111111111111111111111111111111111111111111111111111111111111111111), (83, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111), (89, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111), (97, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111)] 100以下のレプ・ユニット素数の候補が即座に表示できた。 juliaのいいところは,mapやfilterのリストへの操作が簡潔にかけるところだ。 デメリットはPythonのように、何もしなくても多桁(多倍長整数)対応になるわけではなく、 BigIntでくるむ必要があるため、可読性が落ちることだね。 2007年までに9個しか見つかっていないらしい。 pythonで検証してみる。 2,19の2個だけで数分かかる。 次は23の1個だけでもかなり時間がかかる。 317を見つけるには何かしら手段をさらに工夫が必要となりそうだ。 [IN]Python #=========================================================== def isPrime(n): lim = int(round(n**0.5)) for num in range(2,lim): if n % num ==0: return False return True M=list(range(2,20)) divCount=lambda n:len(list(filter(lambda m:n % m==0, list(range(1,n+1))))) P=[x for x in M if divCount(x)==2] R=list(map(lambda n:(n,((10**n-1)//9)),P)) Rp=list(filter(lambda n:isPrime((10**n-1)//9),P)) print(P) print(R) print(Rp) #=========================================================== [2, 3, 5, 7, 11, 13, 17, 19] [(2, 11), (3, 111), (5, 11111), (7, 1111111), (11, 11111111111), (13, 1111111111111), (17, 11111111111111111), (19, 1111111111111111111)] [2, 19] geogebraではn=19けた以上は探さずに終了している。 n=19というので、juliaが19けたが桁溢れしたことを思い出すね。 juliaでも同様なコードをかいて実行しても、時短にはならなかった。残念。 質問:juliaでn=19までのレプ・ユニット素数を検索するにはどうかいたらよいでしょうか。