技術メモ

役に立てる技術的な何か、時々自分用の覚書。幅広く色々なことに興味があります。

Nが現れる素数(N=1,2,3,4)

2が現れる素数という面白い素数が紹介されていた。
2が現れる素数 - INTEGERS

昔せっかく高速素数判定器を作ったので、どうせならNが現れる素数を見つけてやろう!と思い立った。

プログラム

(※プログラムはpython(2.7.12)で動作します)

ルールとしては
①四隅のみの数字を変える(もちろん先頭は1以上の数字)
②四隅の数字はN以外の数字にする
としています。

なので、それぞれ5832(8*9*9*9)個の数字の中から素数を探すことになります。

高速素数判定のプログラム(再掲)

primechecker.pyという名前で保存

import random
import numpy as np

class PrimeChecker:
    def __init__(self, list_limit = pow(10,3)):
        if list_limit < 5: list_limit = 5
        self.prime_list = self.set_prime_list(list_limit)
        self.list_limit = list_limit

    @staticmethod
    def set_prime_list(n):
        def mark(primes,x):
            for i in xrange(2*x, n+1, x):
                primes[i-2] = False
        primes = [(i%2!=0 and i%3!=0 and i%5!=0) for i in range(2,n+1)]
        primes[0] = primes[1] = primes[3] = True
        for x in xrange(7, int(np.sqrt(n+1))+1):
           if primes[x-2]: mark(primes, x)
        prime_list = []
        for i in range(n-1):
            if primes[i]: prime_list.append(i+2)
        return prime_list

    def is_prime(self, n):
        k = 50 # number of iteration times
        if n < self.list_limit: return True if n in self.prime_list else False
        for prime in self.prime_list:
            if n%prime == 0: return False
        else:
            s, d = 0, n-1
            while d%2==0: s,d = s+1, d/2
            while k > 0:
                k = k-1
                a = random.randint(1,n-1)
                t, y = d, pow(a,d,n)
                while t != n-1 and y != 1 and y != n-1:
                    y = pow(y,2,n)
                    t <<= 1
                if y != n-1 and t%2 == 0:
                    return False
            return True

if __name__=="__main__":
    pc = PrimeChecker()
    pc.is_prime(58427) # True

Nが現れる素数を見つけるためのプログラム

from primechecker import PrimeChecker

def display_letter(p):
    """
    display letter on commandline for checking
    """
    print("-"*18)
    for i in range(0,216,18)[::-1]:
        print(str(p/pow(10,i) + pow(10,18))[1:])
        p = p%pow(10,i)
    print("-"*18)

def tex_command(p,num,f):
    """
    output tex command to file
    """
    f.write("\\begin{align} \n")
    for i in range(0,216,18)[::-1]:
        f.write("&")
        string = str(p/pow(10,i) + pow(10,18))[1:]
        string = string.replace(str(num), "\color{red}{%d}"%num)
        f.write(string)
        p = p%pow(10,i)
        f.write("\\\\ \n")
    f.write("\\end{align} \n")

letters = ["""
000000000000000000
000000001111000000
000000011111000000
000000111111000000
000001101111000000
000000001111000000
000000001111000000
000000001111000000
000000001111000000
000000001111000000
000001111111111000
000000000000000000
""","""
000000000000000000
000000222222000000
000022222222220000
000222000002220000
000000000022220000
000000000222200000
000000002222000000
000000022220000000
000002222000000000
000222222222222000
000222222222222000
000000000000000000
""","""
000000000000000000
000000333333000000
000033333333330000
000333000003333000
000000000033333000
000000033333300000
000000033333300000
000000000033333000
000333000003333000
000033333333330000
000000333333000000
000000000000000000
""","""
000000000000000000
000000004444400000
000000044444400000
000000444444400000
000004440444400000
000044400444400000
000444000444400000
004444444444444400
004444444444444400
000000000444400000
000000000444400000
000000000000000000
"""]

pc = PrimeChecker()

f = open("result.txt","w")
for n, letter in enumerate(letters):
    num = n+1
    prime_letter = int(letter.replace("\n",""))
    for i in range(1000,10000)[::-1]:
        if (i/1000)!=num and ((i%1000)/100)!=num and ((i%100)/10)!=num and i%10 != num:
            pl = prime_letter
            pl += pow(10,215) * (i/1000)
            pl += pow(10,198) * ((i%1000)/100)
            pl += pow(10,17) * ((i%100)/10)
            pl += i%10
            if pc.is_prime(pl):
                display_letter(pl)
                tex_command(pl,num,f)
                break
f.close()

実際の素数

1

900000000000000007000000001111000000000000011111000000000000111111000000000001101111000000000000001111000000000000001111000000000000001111000000000000001111000000000000001111000000000001111111111000400000000000000003
\begin{align}
&900000000000000007\\
&00000000\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}000000\\
&0000000\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}000000\\
&000000\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}000000\\
&00000\color{red}{1}\color{red}{1}0\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}000000\\
&00000000\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}000000\\
&00000000\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}000000\\
&00000000\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}000000\\
&00000000\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}000000\\
&00000000\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}000000\\
&00000\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}\color{red}{1}000\\
&400000000000000003\\
\end{align}

2

900000000000000004000000222222000000000022222222220000000222000002220000000000000022220000000000000222200000000000002222000000000000022220000000000002222000000000000222222222222000000222222222222000400000000000000009
\begin{align}
&900000000000000004\\
&000000\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}000000\\
&0000\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}0000\\
&000\color{red}{2}\color{red}{2}\color{red}{2}00000\color{red}{2}\color{red}{2}\color{red}{2}0000\\
&0000000000\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}0000\\
&000000000\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}00000\\
&00000000\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}000000\\
&0000000\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}0000000\\
&00000\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}000000000\\
&000\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}000\\
&000\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}\color{red}{2}000\\
&400000000000000009\\
\end{align}

3

800000000000000002000000333333000000000033333333330000000333000003333000000000000033333000000000033333300000000000033333300000000000000033333000000333000003333000000033333333330000000000333333000000400000000000000009
\begin{align}
&800000000000000002\\
&000000\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}000000\\
&0000\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}0000\\
&000\color{red}{3}\color{red}{3}\color{red}{3}00000\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}000\\
&0000000000\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}000\\
&0000000\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}00000\\
&0000000\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}00000\\
&0000000000\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}000\\
&000\color{red}{3}\color{red}{3}\color{red}{3}00000\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}000\\
&0000\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}0000\\
&000000\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}\color{red}{3}000000\\
&400000000000000009\\
\end{align}

4

900000000000000009000000004444400000000000044444400000000000444444400000000004440444400000000044400444400000000444000444400000004444444444444400004444444444444400000000000444400000000000000444400000500000000000000003
\begin{align}
&900000000000000009\\
&00000000\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}00000\\
&0000000\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}00000\\
&000000\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}00000\\
&00000\color{red}{4}\color{red}{4}\color{red}{4}0\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}00000\\
&0000\color{red}{4}\color{red}{4}\color{red}{4}00\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}00000\\
&000\color{red}{4}\color{red}{4}\color{red}{4}000\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}00000\\
&00\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}00\\
&00\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}00\\
&000000000\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}00000\\
&000000000\color{red}{4}\color{red}{4}\color{red}{4}\color{red}{4}00000\\
&500000000000000003\\
\end{align}

本当は5,6,7,8,9と求めたかったけれど、文字を作るのが面倒で断念してしまいました。
誰か頑張って作ってください。。。

(追記)

そこまで珍しいことだとは思わず、一例しか挙げていませんでした。
存在自体が不思議という話があったので、数を数えてみたのですが
1: 14個
2: 15個
3: 4個
4: 17個
と、そこまで珍しくないことがわかります*1

素数定理に従えば、216桁以下の素数の存在確率は1/log(10^217)≒0.2%らしいです。
リーマン予想の意味,素数分布との関係 | 高校数学の美しい物語
5832*0.002 = 11.6個
だいたいそれくらいですね。

5832個の数字が独立に0.2%の確率で素数かどうか決まるとすると、
各Nに対してそのような素数が存在しない確率は実に0.000850%と非常に低い確率であることがわかります。

ASCIIコード素数を作りたい

それでは、画素数を増やしたいとすればどうなるのでしょうか。
ここでASCIIコードの文字(93文字分)を出来るだけ大きい画素数素数で表現したいとします。
「各文字に対して99.9%の確率でそのような素数が作れます」と言われれば、まあ使い物になりそうだなと言えそうですね。

素数定理に従えば、桁数をnとした時、
\displaystyle
\left( 1 - \frac{1}{nlog10} \right)^{5832} > 0.001
を満たす最大のnを求めればよいことになります。
これを計算するとn=366桁で、
「ASCIIコード素数を作るには19*19くらいの画素数まで大きくできそうだ」
ということが分かりますね。

是非誰か19×19の素数でASCIIコード素数を作ってみてください

*1:4~17/5832を珍しいと考えるかどうかは人によるけれど

2次元配列の特定の条件を満たす要素の中からランダムに選択【numpy】

「0~1の値を取る2次元配列があって、0.7以上の要素を持つインデックスの中からランダムにインデックスを取ってくる。」
こんなことをしたい時にnumpyを駆使すれば綺麗に書けたのでメモしておく。


例えばこんな配列(numpyのarray)があったとする

import numpy as np
a = np.random.randint(0,10,size=[8,6])
array([[6, 9, 4, 1, 2, 4],
       [4, 9, 3, 8, 9, 7],
       [9, 6, 6, 8, 9, 7],
       [6, 6, 2, 7, 2, 3],
       [9, 9, 3, 2, 1, 0],
       [9, 9, 0, 3, 8, 2],
       [0, 5, 4, 4, 0, 7],
       [1, 8, 6, 0, 9, 4]])

この中から8以上の値を持つ要素は

a > 7
array([[False,  True, False, False, False, False],
       [False,  True, False,  True,  True, False],
       [ True, False, False,  True,  True, False],
       [False, False, False, False, False, False],
       [ True,  True, False, False, False, False],
       [ True,  True, False, False,  True, False],
       [False, False, False, False, False, False],
       [False,  True, False, False,  True, False]], dtype=bool)

これだけある。
このTrueの中からランダムに一つインデックスを選びたい。
例えば、[0,1]とかだ。

これをnumpyを使わずに書こうとするとこんな風になるだろう。

import random
index_list = []
for i in range(len(a)):
    for j in range(len(a[0])):
        if a[i][j] > 7:
            index_list.append([i,j])
coord = random.choice(index_list)

randomの便利な関数のおかげでこれでも割とすっきりしているが、
numpyを使えばさらにすっきり書ける。

index = np.random.choice(np.where(a.reshape(-1,) > 7)[0])
coord = [index/a.shape[1], index%a.shape[1]]

たったこれだけだ。
numpy凄い。



※解説

  • reshape[-1,] ... 2次元配列を1次元配列に変換できる
  • np.where ... Trueの箇所のインデックスを返してくれる
  • np.random.choice ... 配列の中から一つ選び取ってくれる

slackbotでbotに投げられた画像をダウンロードする

slackbotでシステムを作ろうとした時に、botに投げた画像をいったん保存する方法がわからなかったので、記録しておく。

slackbot自体の始め方

slackbot自体の使い方は他のブログなどで詳しく紹介されているのでそちらを参考にしてほしい。
インストールはこれだけ

pip install slackbot

自分がslackbotを作る際参考にさせていただいたブログはこちら。
凄く丁寧でわかりやすかった。ありがとうございました。
http://blog.bitmeister.jp/?p=3892

画像の保存

botに画像を投げて色々やらせたかったのだが、どうも保存がうまくいかず手間取ってしまった。
普通に画像をURLからダウンロードする方法ではうまくいかなかった。

import urllib
url = message.body['file']['url_private'] # messageにはファイル付きのメッセージが入っている
urllib.urlretrieve(url)

この方法では数十KBの画像ファイルができるだけで、何も保存できていなかった。

色々調べた結果、GETメソッドを使う時ヘッダーを付けなければいけないらしいとわかった。
(slack公式サイト:file type | Slack)
正直この辺はよくわかっていないんだけど、requestsを使えば何とかなるらしい。

ということでサンプルがこちら。
プラグインのファイルの中にtokenをオリジナルのものにしてこれを書く。

import requests

@listen_to("(.*)")
def img(message, params):
    if 'file' in message.body:
        url = message.body['file']['url_private']
        flag = message.body['file']['filetype']
        tmpfile = "./tmp." + flag

        token = 'YOUR_TOKEN'
        rst = requests.get(url, headers={'Authorization': 'Bearer %s' % token}, stream=True)

        fo = open(tmpfile, "wb")
        shutil.copyfileobj(rst.raw, fo)
        fo.close()


今回はpythonのslackbotを紹介したけれど、この問題だけなら他のライブラリや言語でもこれを応用して解決できるかもしれない。