Pebble Coding

プログラマーの作業メモ

iOS11 Audio関連機能(2017)

iOS11 Audio関連で追加された機能で気になったものをまとめておきます。
内容はここで見ることができます。
今年からApple Developer登録していない方でもVideoが見られるようになったようです。
(WWDC 2017に見た新しい教材のあり方|MacFan)

developer.apple.com

AVAudioEngine関連

offlineモード追加

ようやくofflineモードが追加され、オーディオファイルに出力できるようになりました。
Mixdownとも言います。
CoreAudioでは元々できましたが、AVAudioEngineでもできるようになりました。

realtimeレンダリング

サンプル単位でリアルタイムでデータを入れて再生できるようになりました。
CoreAudioで実装する場合と同じくUnsafePointer<AudioBufferList>を使うようです。

[補足]動作するサンプルコードがどこにも見当たりません。

Completion Callbacks

.dataConsumed
.dataRendered
.dataPlayedBack

いまいち 理解できておりませんがメモ。

AUGraphが2018年にDeprecationとなる予定

ううむ、マジですか。
iOS12にはAUGraphの全ての機能をAVAudioEngineに入れるつもりってことですな。
まあ、Objective-Cでそのまま書き換えればいいかな。
すぐに使えなくなるわけでもないので放置でもいいですが。

AirPlay2関連

AirPlay自体よく分かっておりません。。

AUAudioUnit関連

MIDI出力できるようになったようです。
iOSにパラメータ表示Viewが追加されたようです。

FLAC,Opus に対応

FLACは可逆フォーマット(wav,cafと同じ)
Opusは非可逆フォーマット(aac,mp3と同じ)

ハイレゾプレーヤーを作っている方は
自前で実装しなくてよくなったということですな。

Iner-Device Audio Mode + MIDI

macOS 10.11(El Capitan),iOS9でAudioのやりとりができるように追加された機能ですが、
macOS 10.13(High Sierra),iOS11でMIDIのやりとりもできるようになったようです。

デモでは、iPadシーケンサーからMIDI(ノートやボリュームオートメーション?)
を出力しmacOSのLogicで録音していました。

指定の長さのsin波のWAVファイルの作り方(モノラル、ステレオ)

指定の長さのsin波のWAVファイルを簡単に作る方法をメモしておきます。

Audacityというフリーのアプリケーションを使います。

Audacity® | Free, open source, cross-platform audio software for multi-track recording and editing.

Mac/Windows/Linux版があります。

Mac版は画面のここからダウンロードします。

f:id:pebble8888:20170914195711p:plain:w600

Audacityをインストールして立ち上げるとこのような画面になります。

f:id:pebble8888:20170914200122p:plain:w600

Sin波を作るには [ジェネレータ] - [トーン] を選びます。

f:id:pebble8888:20170914200422p:plain:w200

ここでは、
波形: サイン波
周波数: 440
振幅: 1
継続時間: 10s

を選びました。

f:id:pebble8888:20170914200555p:plain:w400

10秒のモノラルの音声データが生成されました。
再生ボタンを押して再生して確かめることができます。

f:id:pebble8888:20170914200757p:plain:w600

ファイルに保存するには、[ファイル]-[オーディオの書き出し]を選択して、
ファイルタイプは[WAV (Microsoft) 16 bit PCM 符号あり]を選び保存したら完了です。

ステレオファイルの作成

モノラルファイルの作り方は分かったので、次はステレオファイルを作ってみます。
先ほどモノラルトラックを作った状態からさらに、[トラック]-[新しく追加]-[モノラルトラック]を選択します。

f:id:pebble8888:20170914202310p:plain:w600

同じようにジェネレーターを使ってサイン波を生成します。

f:id:pebble8888:20170914202327p:plain:w600

作成した2つのモノラルトラックを一つにまとめるため、[編集]-[選択]-[全て]で選択状態にします。

f:id:pebble8888:20170914202347p:plain:w600

次に分かりづらいですが、ここの三角形をクリックします。

f:id:pebble8888:20170914202619p:plain:w200

[ステレオトラックの作成]を選択します。
これであとはファイルに保存するだけです。

f:id:pebble8888:20170914202634p:plain:w200

swift であらゆるデータをHex表示してデバッグする

Xcodeのデバッガではswiftオブジェクトの中身を表示しても生のデータを表示してくれず、デバッグしづらいです。
(C/C++なら大丈夫ですが。)
あらゆるものを16進Hex表示の文字列にしてデバッグに使いましょう。

public protocol HexRepresentable {
    func hexDescription() -> String
}

こんなプロトコルでいいでしょう。

extension UInt8 : HexRepresentable
{
    public func hexDescription() -> String {
        return String(format:"%02x", self)
    }
}

extension Int32 : HexRepresentable
{
    public func hexDescription() -> String
    {
        return String(format:"%08x", self)
    }
}

extension Collection where Iterator.Element : HexRepresentable {
    public func hexDescription() -> String {
        return self.map({ $0.hexDescription() }).joined()
    }
}

extension Data : HexRepresentable {
    public func hexDescription() -> String {
        return self.map({ String(format:"%02x", $0) }).joined() 
    }
}

これだけあれば、事足りるでしょう。

参照:

PBKit/PBHex.swift at master · pebble8888/PBKit · GitHub

楕円曲線上の離散対数問題(ECDLP)

一つ前の記事

python で有限体Fpでの楕円曲線上の有理点の群構造を調べる - Pebble Coding

で有限体F_{p}での楕円曲線の群構造をpythonを使って実験的に調べました。
分かったことは、
- 群の位数=群の元の個数=有理点の個数(無限遠点を含む)である。
- 群の位数は必ずしも素数ではない、むしろそうでない方が多い。
- 有理点の位数は群の位数の約数になるっぽい。

ここで、
有理点の位数(もしくは群の元の位数)の定義は、
有理点Pをn個加算した時に初めてO(無限遠点)になることです。

楕円曲線上の離散対数問題(Elliptic Curve Discrete Logarithm Problem)を定義します。

離散対数問題の解説はこちら 離散対数問題(Discrete Logarithm Problem)とは - Pebble Coding

有限体F_{p}上の楕円曲線の2つの元(有理点)であるH,Gに対して、
H=xG (つまりGをx回加算する)となるxが存在するならxを求めよ。
という問題です。

各種攻撃法はいくつかありますが、攻撃を回避するように注意して設定すれば、
問題ありません。

この群の位数(有理点の数)は必ずしも素数である必要はありません。
群の位数Nが例えばN=8rだとします。
rは2^{160}くらいの大きな素数とします。

有理点の位数は群位数の約数になるので、
2, 4, 8, 2r, 4r, 8r
の6パターンのどれかとなりますが、2, 4, 8 の点を選ぶのは論外なので、
2r, 4r, 8r のどれかの点をGとして選びます。
そうすれば、ECDLPの強度としては十分だと言えます。

実際に
Edwards-Curve Digital Signature Algorithm (EdDSA)
のRFC8032ではベースポイントGをこのように選ぶようにと記載されています。
RFC 8032 - Edwards-Curve Digital Signature Algorithm (EdDSA)

このような楕円曲線をどのように選べば良いのでしょうか?
加算の計算が高速にできるものを選ぶ必要もあります。
この問題は暗号学者でない我々にとってはあまり気にする必要はありません。
というのも、現在普及しているEd25519は、
楕円曲線の形、有限体、ベースポイント(ベースポイントの位数)、の全てがほぼ決まっておりRFCとなっているからです。
opensslなどにはすでに実装されていますし、各種言語のライブラリにもそのうち実装されるでしょう。

参考:

代数学から学ぶ暗号理論: 整数論の基礎から楕円曲線暗号の実装まで

代数学から学ぶ暗号理論: 整数論の基礎から楕円曲線暗号の実装まで

暗号理論入門 原書第3版

暗号理論入門 原書第3版

楕円曲線論入門

楕円曲線論入門

python で有限体Fpでの楕円曲線上の有理点の群構造を調べる

ここでは、有限体F_{5}
(p=5)
楕円曲線 y^{2} = x^{3} + x + 1
(a=0,b=1,c=1) の有理点をpythonで調べています。 有理点の数は9です。(無限遠点を含む)
無限遠点はOと出力しています。
加法公式を用いて、有理点{P1, P1, … P8}を2倍,3倍,…,9倍した点も示しています。
この計算の途中で元の有理点に戻った場合、あとは繰り返しなので、そこで計算を止めています。

P1,P2,P5,P6,P7,P8は9倍したところで無限遠点になるので、位数9の点といいます。 ここで倍々していくと9つ全ての点を通ります。
これを巡回群であるといいます。

P3,P4は3倍したところで無限遠点になるので、位数3の点といいます。

この群は2つの巡回群Z/3ZとZ/9Zの直積と同型であるといいます。

%matplotlib inline 
import numpy as np
import matplotlib.pyplot as plt

class Point:
    def __init__(self, x, y):
        self.x = x
        self.y = y
    
    def desc(self):
        if self.x == -1 and self.y == -1:
            return "O"
        return "("+str(self.x)+","+str(self.y)+")"
    
    def is_equal(self, other):
        return self.x == other.x and self.y == other.y

    def is_equal_negative(self, other, prime):
        return self.x == other.x and - self.y == other.y - prime
        
    def iszero(self):
        return self.x == -1 and self.y == -1
    
    @staticmethod
    def zero():
        return Point(-1, -1)

class EC:
    def __init__(self, a, b, c, p):
        self.a = a
        self.b = b
        self.c = c
        self.p = p
        self.calc_order()
    
    def calc_order(self):
        self.points = []
        for x in range(self.p):
            for y in range(self.p):
                pt = Point(x,y)
                if self.oncurve(pt):
                    self.points.append(pt)
        self.points_count = len(self.points)
        self.order = len(self.points)+1

    def oncurve(self, pt):
        l = ((pt.y**2) % self.p)
        r = ((pt.x**3) + self.a*(pt.x**2) + self.b*pt.x + self.c) % self.p
        return l == r

    def d(self):
        return -4*(self.a**3)*self.c \
                + (self.a**2)*(self.b**2) \
                + 18*self.a*self.b*self.c \
                - 4*(self.b**3) \
                - 27 * (self.c**2)
    
    def inverse(self, a):
        assert self.p >= 2
        return a ** (self.p-2)
    
    def plus(self, p1, p2):
        x1 = p1.x
        y1 = p1.y
        x2 = p2.x
        y2 = p2.y
        if (p1.iszero()):
            return p2
        if (p2.iszero()):
            return p1
        
        if p1.is_equal_negative(p2, self.p):
            return Point.zero()
        
        elif p1.is_equal(p2):
            if (y1 == 0):
                return Point.zero()
            else:
                lm = (3 * (x1 **2) + 2 * self.a * x1 + self.b) * self.inverse(2 * y1)
                x3 = (lm ** 2) - self.a - x1 - x2
                y3 = - lm *(x3-x1) - y1
        else:
            lm = (y2-y1) * self.inverse(x2-x1)
            x3 = (lm ** 2) - self.a - x1 - x2
            y3 = -(lm*(x3-x1) + y1)

        q = Point(x3 , y3)
        q2 = Point(x3 % self.p, y3 % self.p)
        if not self.oncurve(q2):
            print("p1 {0} p2 {1} q {2} q2 {3}".format(p1.desc(), p2.desc(), q.desc(), q2.desc()))
            assert False
        return q2

    def mul(self, p, n):
        q = p
        for i in range(n-1):
            q = ec.plus(q, p)
        return q

ec = EC(0, 1, 1, 5)

print "F" + str(ec.p)
print "d:" + str(ec.d())
print "#EC:" + str(ec.order)
print ""

plotx = [p.x for p in ec.points]
ploty = [p.y for p in ec.points]
n = ["P"+str(i+1) for (i, p) in zip(range(len(ec.points)), ec.points)]

for i in range(ec.points_count):
    baseP = ec.points[i]
    print "P"+str(i+1)+":" + baseP.desc()
    for j in range(2, ec.order+1):
        p = ec.mul(baseP, j)
        print str(j)+"*P"+str(i+1)+":"+p.desc()
        if (p == baseP):
            print ""
            break
    print ""
    
fig, ax = plt.subplots()
ax.scatter(plotx, ploty)

for i, txt in enumerate(n):
    ax.annotate(txt, (plotx[i],ploty[i]))
    
F5
d:-31
#EC:9

P1:(0,1)
2*P1:(4,2)
3*P1:(2,1)
4*P1:(3,4)
5*P1:(3,1)
6*P1:(2,4)
7*P1:(4,3)
8*P1:(0,4)
9*P1:O

P2:(0,4)
2*P2:(4,3)
3*P2:(2,4)
4*P2:(3,1)
5*P2:(3,4)
6*P2:(2,1)
7*P2:(4,2)
8*P2:(0,1)
9*P2:O

P3:(2,1)
2*P3:(2,4)
3*P3:O
4*P3:(2,1)


P4:(2,4)
2*P4:(2,1)
3*P4:O
4*P4:(2,4)


P5:(3,1)
2*P5:(0,1)
3*P5:(2,4)
4*P5:(4,2)
5*P5:(4,3)
6*P5:(2,1)
7*P5:(0,4)
8*P5:(3,4)
9*P5:O

P6:(3,4)
2*P6:(0,4)
3*P6:(2,1)
4*P6:(4,3)
5*P6:(4,2)
6*P6:(2,4)
7*P6:(0,1)
8*P6:(3,1)
9*P6:O

P7:(4,2)
2*P7:(3,4)
3*P7:(2,4)
4*P7:(0,4)
5*P7:(0,1)
6*P7:(2,1)
7*P7:(3,1)
8*P7:(4,3)
9*P7:O

P8:(4,3)
2*P8:(3,1)
3*P8:(2,1)
4*P8:(0,1)
5*P8:(0,4)
6*P8:(2,4)
7*P8:(3,4)
8*P8:(4,2)
9*P8:O

f:id:pebble8888:20170909165102p:plain

他の楕円曲線でも試してみましょう。

F_{17}
y^{2}=x^{3}+2 x+4

この場合、全ての点の位数が17となります。こういうのはわりと珍しいです。 この群はZ/17Zと同型であるといいます。

F13
d:-464
#EC:17

P1:(0,2)
2*P1:(10,6)
3*P1:(12,1)
4*P1:(2,9)
5*P1:(7,6)
6*P1:(5,10)
7*P1:(9,7)
8*P1:(8,8)
9*P1:(8,5)
10*P1:(9,6)
11*P1:(5,3)
12*P1:(7,7)
13*P1:(2,4)
14*P1:(12,12)
15*P1:(10,7)
16*P1:(0,11)
17*P1:O

P2:(0,11)
2*P2:(10,7)
3*P2:(12,12)
4*P2:(2,4)
5*P2:(7,7)
6*P2:(5,3)
7*P2:(9,6)
8*P2:(8,5)
9*P2:(8,8)
10*P2:(9,7)
11*P2:(5,10)
12*P2:(7,6)
13*P2:(2,9)
14*P2:(12,1)
15*P2:(10,6)
16*P2:(0,2)
17*P2:O

P3:(2,4)
2*P3:(8,5)
3*P3:(7,6)
4*P3:(0,2)
5*P3:(12,12)
6*P3:(9,6)
7*P3:(5,10)
8*P3:(10,6)
9*P3:(10,7)
10*P3:(5,3)
11*P3:(9,7)
12*P3:(12,1)
13*P3:(0,11)
14*P3:(7,7)
15*P3:(8,8)
16*P3:(2,9)
17*P3:O

P4:(2,9)
2*P4:(8,8)
3*P4:(7,7)
4*P4:(0,11)
5*P4:(12,1)
6*P4:(9,7)
7*P4:(5,3)
8*P4:(10,7)
9*P4:(10,6)
10*P4:(5,10)
11*P4:(9,6)
12*P4:(12,12)
13*P4:(0,2)
14*P4:(7,6)
15*P4:(8,5)
16*P4:(2,4)
17*P4:O

P5:(5,3)
2*P5:(7,6)
3*P5:(0,11)
4*P5:(9,6)
5*P5:(2,9)
6*P5:(10,7)
7*P5:(8,5)
8*P5:(12,1)
9*P5:(12,12)
10*P5:(8,8)
11*P5:(10,6)
12*P5:(2,4)
13*P5:(9,7)
14*P5:(0,2)
15*P5:(7,7)
16*P5:(5,10)
17*P5:O

P6:(5,10)
2*P6:(7,7)
3*P6:(0,2)
4*P6:(9,7)
5*P6:(2,4)
6*P6:(10,6)
7*P6:(8,8)
8*P6:(12,12)
9*P6:(12,1)
10*P6:(8,5)
11*P6:(10,7)
12*P6:(2,9)
13*P6:(9,6)
14*P6:(0,11)
15*P6:(7,6)
16*P6:(5,3)
17*P6:O

P7:(7,6)
2*P7:(9,6)
3*P7:(10,7)
4*P7:(12,1)
5*P7:(8,8)
6*P7:(2,4)
7*P7:(0,2)
8*P7:(5,10)
9*P7:(5,3)
10*P7:(0,11)
11*P7:(2,9)
12*P7:(8,5)
13*P7:(12,12)
14*P7:(10,6)
15*P7:(9,7)
16*P7:(7,7)
17*P7:O

P8:(7,7)
2*P8:(9,7)
3*P8:(10,6)
4*P8:(12,12)
5*P8:(8,5)
6*P8:(2,9)
7*P8:(0,11)
8*P8:(5,3)
9*P8:(5,10)
10*P8:(0,2)
11*P8:(2,4)
12*P8:(8,8)
13*P8:(12,1)
14*P8:(10,7)
15*P8:(9,6)
16*P8:(7,6)
17*P8:O

P9:(8,5)
2*P9:(0,2)
3*P9:(9,6)
4*P9:(10,6)
5*P9:(5,3)
6*P9:(12,1)
7*P9:(7,7)
8*P9:(2,9)
9*P9:(2,4)
10*P9:(7,6)
11*P9:(12,12)
12*P9:(5,10)
13*P9:(10,7)
14*P9:(9,7)
15*P9:(0,11)
16*P9:(8,8)
17*P9:O

P10:(8,8)
2*P10:(0,11)
3*P10:(9,7)
4*P10:(10,7)
5*P10:(5,10)
6*P10:(12,12)
7*P10:(7,6)
8*P10:(2,4)
9*P10:(2,9)
10*P10:(7,7)
11*P10:(12,1)
12*P10:(5,3)
13*P10:(10,6)
14*P10:(9,6)
15*P10:(0,2)
16*P10:(8,5)
17*P10:O

P11:(9,6)
2*P11:(12,1)
3*P11:(2,4)
4*P11:(5,10)
5*P11:(0,11)
6*P11:(8,5)
7*P11:(10,6)
8*P11:(7,7)
9*P11:(7,6)
10*P11:(10,7)
11*P11:(8,8)
12*P11:(0,2)
13*P11:(5,3)
14*P11:(2,9)
15*P11:(12,12)
16*P11:(9,7)
17*P11:O

P12:(9,7)
2*P12:(12,12)
3*P12:(2,9)
4*P12:(5,3)
5*P12:(0,2)
6*P12:(8,8)
7*P12:(10,7)
8*P12:(7,6)
9*P12:(7,7)
10*P12:(10,6)
11*P12:(8,5)
12*P12:(0,11)
13*P12:(5,10)
14*P12:(2,4)
15*P12:(12,1)
16*P12:(9,6)
17*P12:O

P13:(10,6)
2*P13:(2,9)
3*P13:(5,10)
4*P13:(8,8)
5*P13:(9,6)
6*P13:(7,7)
7*P13:(12,12)
8*P13:(0,11)
9*P13:(0,2)
10*P13:(12,1)
11*P13:(7,6)
12*P13:(9,7)
13*P13:(8,5)
14*P13:(5,3)
15*P13:(2,4)
16*P13:(10,7)
17*P13:O

P14:(10,7)
2*P14:(2,4)
3*P14:(5,3)
4*P14:(8,5)
5*P14:(9,7)
6*P14:(7,6)
7*P14:(12,1)
8*P14:(0,2)
9*P14:(0,11)
10*P14:(12,12)
11*P14:(7,7)
12*P14:(9,6)
13*P14:(8,8)
14*P14:(5,10)
15*P14:(2,9)
16*P14:(10,6)
17*P14:O

P15:(12,1)
2*P15:(5,10)
3*P15:(8,5)
4*P15:(7,7)
5*P15:(10,7)
6*P15:(0,2)
7*P15:(2,9)
8*P15:(9,7)
9*P15:(9,6)
10*P15:(2,4)
11*P15:(0,11)
12*P15:(10,6)
13*P15:(7,6)
14*P15:(8,8)
15*P15:(5,3)
16*P15:(12,12)
17*P15:O

P16:(12,12)
2*P16:(5,3)
3*P16:(8,8)
4*P16:(7,6)
5*P16:(10,6)
6*P16:(0,11)
7*P16:(2,4)
8*P16:(9,6)
9*P16:(9,7)
10*P16:(2,9)
11*P16:(0,2)
12*P16:(10,7)
13*P16:(7,7)
14*P16:(8,5)
15*P16:(5,10)
16*P16:(12,1)
17*P16:O

もう一つ試します。無駄に長くなってきたので、点の位数(Pをn倍した時に無限遠点Oになるもっとも小さな値)だけ表示します。
 F_{71}
 y^{2} = x^{3} - x

有理点の数(=群の位数)は72, 有理点の位数は、2, 3, 4, 9, 12, 18, 36で有理点の数の因数が全て出てきているようです。

F71
d:4
#EC:72

P1:(0,0)
order:2


P2:(1,0)
order:2


P3:(2,19)
order:9


P4:(2,52)
order:9


P5:(3,33)
order:18


P6:(3,38)
order:18


P7:(4,29)
order:9


P8:(4,42)
order:9


P9:(5,7)
order:18


P10:(5,64)
order:18


P11:(9,9)
order:6


P12:(9,62)
order:6


P13:(12,15)
order:36


P14:(12,56)
order:36


P15:(13,14)
order:4


P16:(13,57)
order:4


P17:(14,23)
order:18


P18:(14,48)
order:18


P19:(19,33)
order:3


P20:(19,38)
order:3


P21:(21,9)
order:36


P22:(21,62)
order:36


P23:(23,28)
order:18


P24:(23,43)
order:18


P25:(27,29)
order:36


P26:(27,42)
order:36


P27:(32,17)
order:12


P28:(32,54)
order:12


P29:(33,7)
order:36


P30:(33,64)
order:36


P31:(35,13)
order:18


P32:(35,58)
order:18


P33:(37,8)
order:9


P34:(37,63)
order:9


P35:(40,29)
order:12


P36:(40,42)
order:12


P37:(41,9)
order:36


P38:(41,62)
order:36


P39:(42,8)
order:18


P40:(42,63)
order:18


P41:(43,21)
order:36


P42:(43,50)
order:36


P43:(45,22)
order:36


P44:(45,49)
order:36


P45:(46,34)
order:36


P46:(46,37)
order:36


P47:(47,20)
order:18


P48:(47,51)
order:18


P49:(49,33)
order:18


P50:(49,38)
order:18


P51:(51,16)
order:12


P52:(51,55)
order:12


P53:(53,24)
order:18


P54:(53,47)
order:18


P55:(54,28)
order:36


P56:(54,43)
order:36


P57:(55,31)
order:12


P58:(55,40)
order:12


P59:(56,30)
order:6


P60:(56,41)
order:6


P61:(60,10)
order:4


P62:(60,61)
order:4


P63:(61,2)
order:36


P64:(61,69)
order:36


P65:(63,8)
order:6


P66:(63,63)
order:6


P67:(64,27)
order:36


P68:(64,44)
order:36


P69:(65,28)
order:36


P70:(65,43)
order:36


P71:(70,0)
order:2

f:id:pebble8888:20170909172711p:plain

macOS 10.12 jupyer notebook でグラフを描画する(Fpにおける楕円曲線の解の個数)

matplotlibをインストールします。

~$ pip install matplotlib
%matplotlib inline
 
import numpy as np
import matplotlib.pyplot as plt
 
# 乱数を生成
x = np.random.rand(100)
y = np.random.rand(100)
 
# 散布図を描画
plt.scatter(x, y)
plt.show()

うまくいきました。

f:id:pebble8888:20170908220013p:plain

次に y = x^{2}となる整数を少しプロットしてみます。

%matplotlib inline 
import numpy as np
import matplotlib.pyplot as plt


a = np.array([[1,1], [2,4], [3,9], [4, 16]])
x = a[:, 0]
y = a[:, 1]
print x
print y

plt.scatter(x, y)
plt.show()

f:id:pebble8888:20170908222813p:plain

今度は数式で計算します。

%matplotlib inline 
import numpy as np
import matplotlib.pyplot as plt

a = np.empty((0,2))

def square(x):
    return x*x
    
for i in range(5):
    a = np.append(a, np.array([[i,square(i)]]), axis=0)

x = a[:, 0]
y = a[:, 1]
plt.scatter(x, y)
plt.show()

f:id:pebble8888:20170908230438p:plain

さてついに、有限体 F_p上の楕円曲線の解をプロットしてみましょう。
 F_5,  y^{2} = x^{3} + x^{2} + 1の場合です。

%matplotlib inline 
import numpy as np
import matplotlib.pyplot as plt

points = np.empty((0,2))

class EC:
    def __init__(self, a, b, c, p):
        self.a = a
        self.b = b
        self.c = c
        self.p = p
    
    def oncurve(self, x, y):
        l = ((y**2) % self.p)
        r = (((x**3) + self.a*(x**2) + self.b*x + self.c) % self.p)
        #print "(" + str(x) + "," + str(y) + ") " + "l=" + str(l) + ",r=" + str(r)
        return l == r

    def d(self):
        return -4*(self.a**3)*self.c \
                + (self.a**2)*(self.b**2) \
                + 18*self.a*self.b*self.c \
                - 4*(self.b**3) \
                - 27 * (self.c**2)

ec = EC(0, 1, 1, 5)

for x in range(ec.p):
    for y in range(ec.p):
        if ec.oncurve(x, y):
            points = np.append(points, np.array([[x,y]]), axis=0)

print "d:" + str(ec.d())
print "#EC:" + str(points.shape[0]+1)

plotx = points[:, 0]
ploty = points[:, 1]
plt.scatter(plotx, ploty)
plt.show()

f:id:pebble8888:20170908234804p:plain

有理点の数が#EC: 8+1(無限遠点)=9 であることが分かりました。

これで、 F_{p}における任意の楕円曲線上の有理点を簡単に数えられる実験道具を手に入れました。

Hasse-Weilの定理が正しい事を、この曲線  y^{2} = x^{3} + x + 1 で試してみましょう。

f:id:pebble8888:20170909000043p:plain
f:id:pebble8888:20170909000052p:plain
f:id:pebble8888:20170909000945p:plain
f:id:pebble8888:20170909001055p:plain
f:id:pebble8888:20170909001153p:plain
f:id:pebble8888:20170909001257p:plain

確かにHasse-Heilの定理が成り立っている気がしますね。

B-smooth と B-power-smooth

B-smooth と B-power-smooth の概念を押さえておきます。

B-smooth

任意の自然数Nを素因数に分解すると次の形になります。
 N = p_{1}^{e_1} p_{2}^{e_2} p_{3}^{e_3}  \cdot \cdot \cdot p_{i}^{e_i}  \cdot \cdot \cdot p_{m}^{e_m}
ここで、 p_{k}はk番目の素数 {2, 3, 5, 7, 11, …}を表します。
 e_{k}は0以上の整数です。
ここで最も大きい素数 p_m =Bとして、
NをB-smoothであるといいます。

例:
16=2^{4}は2-smoothです。
 9800 = 2^{3} 5^{2} 7^{2}は7-smoothです。

B-power-smooth

任意の自然数Nを素因数に分解すると次の形になります。
 N = p_{1}^{e_1} p_{2}^{e_2} p_{3}^{e_3}  \cdot \cdot \cdot p_{i}^{e_i}  \cdot \cdot \cdot p_{m}^{e_m}
ここで、 p_{k}はk番目の素数 {2, 3, 5, 7, 11, …}を表します。
 e_{k}は0以上の整数です。
このm個の値  p_{1}^{e_1}, p_{2}^{e_2},  p_{3}^{e_3},  \cdot \cdot \cdot,  p_{i}^{e_i},   \cdot \cdot \cdot , p_{m}^{e_m} のうちの最大の値をBとする時、
NをB-power-smoothであるといいます。

例:
 2520= 2^{4} 3^{2} 5 \cdot 7で、  2^{4}, 3^{2}, 5, 7のうち最も大きな値は 2^{4} = 16ですから、
2520は16-power-smoothであるといいます。

NがB-power-smoothであるとき、NはBで割り切れることを覚えておいてください。

参考:

Smooth Number -- from Wolfram MathWorld

Smooth number - Wikipedia