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

[b][u][size=150][color=#9900ff]質問：rustでn=317までのレプ・ユニット素数を検索するには、どうしましょう。[br][/color][br]素数リスト、レピュニット数リスト、レピュニット素数リストの順に作るのは同じですが、[br]最後の素数判定でスピードアップの工夫をします。[br]rustが高速なのと、素数判定が高速なので、あっという間に結果がでます。[br][br][/size][/u][/b]約数2個の整数、つまり、素数のリストをつくりましょう。[br]2からlast＝330までの範囲を対象に[br]約数2個というフィルターをかけて集めます。[br][br] [color=#0000ff] let prime: Vec = (2..=last)[br] .filter(|&n| (1..=n).filter(|&m| n % m == 0).count() == 2)[br] .collect();[br] [br] println!("素数リスト: {:?}", prime);[br][/color][br]次に１が素数個ならんだレピュニット数のリストを作りましょう。[br]素数を対象に、10の素数乗ー１を9で割った商に変換して、集めます。[br] [color=#0000ff] // Rn = (10^n - 1) / 9[br] let repunit: Vec = prime[br] .iter()[br] .map(|&x| (ten.pow(x as u32) - &one) / &nine)[br] .collect();[br][br] println!("レピュニット数リスト（Rn）: {:?}", repunit);[br][/color][br]最後に、レピュニット数Rnが素数かどうかを判定しますが。[br]１が330個ならんだ数が素数かどうかまで調べる計算は、[br][b]通常の素数判定[/b]（約数２個とか、２から平方根までに約数が1個もないとか）では[br]時間が膨大になり、[b]PCが永遠の眠り[/b]につきます。[br]そこで、[b]乱択法[/b]を使いましょう。[br]円の面積をモンテカルロ法で調べるように、乱数を使って近似的に調べます。[br][b]Miller-Rabin法[/b]です。[br][br][b]nが素数ならば、[br][/b]n - 1を２で割り切れるまで何度も割りましょう。[br]その割り算回数をk、最後の奇数の商をqとすると、n-１ = 2^k * q です。[br]2 <= a < n-1 の範囲でランダムな整数 a を選択して、[br][b]a^q ≡ 1 (mod n) または、[br]a^[(2^(i-1))*q] ≡-1(mod n)となるiが１からkの中に１つある。[br][/b][br]このことの対偶から、[br][br]ランダムな整数 aで、mod nで [br][b]a^q !≡ 1 (mod n) and 1からkのどのiでもa^[(2^(i-1))*q] ≡-1(mod n)[br][/b]ならば、[b]nが合成数。このときのaを合成数の証人[/b]といいます。[br][br][b]Miller-Rabin法[/b]では証人aが4分3の合成数以上確保できるので、10回連続失敗する可能性は100万分の１程度と少ない（しかし、ゼロではない。）[br][b][color=#0000ff]合成数の判定がでないと、相当の確率で素数と言える。しかし、素数でない可能性もの残る。[br][/color][/b][br]アルゴリズムとしては、≡1, ≡-1のすべての判定をすべて不合格なら[br]最後に合成数と言えるので、「たぶん素数か」の問いに対してfalseを返す。[br]反対に、≡1, ≡-1が言えた時点で次の試行に進みどれも合格したときにtrueにたどり着く。[br][br][color=#0000ff][b][IN]rust in Cargo.toml [br][dependencies][br]rand = "0.8"[br]num-bigint = { version = "0.4", features = ["rand"] } # rand機能を有効にする[br]num-traits = "0.2"[br]num-iter = "0.1"[br][br][/b][IN]rust in main.rs[br]use num_bigint::{BigInt, RandBigInt, ToBigInt};[br]use num_traits::{One, Zero};[br]use rand::thread_rng; // rng ではなく thread_rng[br][br]fn [b]main[/b]() {[br]    let one: BigInt = One::one();[br]    let ten: BigInt = 10.to_bigint().unwrap();[br]    let nine: BigInt = 9.to_bigint().unwrap();[br][br]    let last = 330;  //n=330までやってみる。[br]    let [b]prime[/b]: Vec< i64 > = (2..=last)[br]        .filter(|&n| (1..=n).filter(|&m| n % m == 0).count() == 2)[br]        .collect();[br]    [br]    println!("素数リスト: {:?}", prime);[br][br]    // Rn = (10^n - 1) / 9[br]    let [b]repunit[/b]: Vec< BigInt > = prime[br]        .iter()[br]        .map(|&x| (ten.pow(x as u32) - &one) / &nine)[br]        .collect();[br][br]    println!("レピュニット数リスト（Rn）: {:?}", repunit);[br][br]    let [b]repunit_prime[/b]: Vec< BigInt > = repunit[br]        .into_iter()[br]        .filter(|rn| is_prime_miller_rabin(rn, 10)) // 10回テスト[br]        .collect();[br][br]    println!("レピュニット素数: {:?}", repunit_prime);[br]}[br][br]fn [b]is_prime_miller_rabin[/b](n: &BigInt, repeat: usize) -> bool {[br]    let zero = BigInt::zero();[br]    let one = BigInt::one();[br]    let two = 2.to_bigint().unwrap();[br][br]    if n <= &one { return false; }[br]    if n <= &3.to_bigint().unwrap() { return true; }[br]    if n % &two == zero { return false; }[br][br]    // n - 1 = 2^k * q を計算して、kとqを絞り出す。[br]    let n_minus_one = n - &one;[br]    let mut q = n_minus_one.clone();[br]    let mut k = 0u32;[br]    while &q % &two == zero {[br]        q /= &two;[br]        k += 1;[br]    }[br]    //乱択の準備[br]    let mut rng = thread_rng(); [br]    //repeat回試行する。[br]    'outer: for _ in 0..repeat {[br]        // 2 <= a < n-1 の範囲でランダムな整数 a を選択[br]        let a = rng.gen_bigint_range(&two, &n_minus_one);[br]        [br]        // v = a^q mod n　を求めてテストの準備[br]        let mut v = a.modpow(&q, n);[br][br]        // a^q ≡ 1 (mod n) なら次の試行へ[br]        if v == one || v == n_minus_one {[br]            continue 'outer;[br]        }[br][br]        // a^{ (2^i) * q } ≡ -1 (mod n) になる i があるかチェック[br]        for _ in 0..k - 1 {[br]            v = v.modpow(&two, n); // 順次 2乗していく[br]            if v == n_minus_one {[br]                continue 'outer;[br]            }[br]        }[br]        // どの条件も満たさなければ、間違いなくnは合成数であり、aはその証人となる。[br]        return false;[br]    }[br]    //たぶん素数[br]    true[br]}[br][/color][OUT]cargo run --quiet[br]素数リスト: [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, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317][br]レピュニット数リスト（Rn）: [11, 111, 11111, 1111111, 11111111111, 1111111111111, 11111111111111111, 1111111111111111111, 11111111111111111111111, 11111111111111111111111111111, 1111111111111111111111111111111, 1111111111111111111111111111111111111, 11111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111][br][b][color=#0000ff]レピュニット素数: [11, 1111111111111111111, 11111111111111111111111, 11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111][br][/color][/b][br][br][br]