🔖

ランダムに3点とって鋭角三角形となる確率をJulia で考える。

2021/07/15に公開

Julia

Twitter

math

idea

はじめに

まずは，次のようなツイートがありました。

いろいろ考えてみることにしました。

確率0

まずは，全平面でランダムに3点取る場合です。

任意に2点を取ると，3点目をどこに取るかということになります。ここで3点目を落とす場所によって，その面積比で確率を考えることにしました。

両方とも面積は $\infty$ ですが，鋭角三角形を作る部分の方が小さいと考えられます。

確率は0としました。

「領域をある程度限定した方がいいのかな？」と思ったので，正方形と円で考えることにしました。

1x1の正方形の内部に3点を任意に取る。

続いて，1x1の正方形の内部に3点を取ります。まずは，100万回作ってその割合を調べることにしました。

square.jl

n=1000000
k=0
for i=1:n
    a1=rand(2)
    a2=rand(2)
    a3=rand(2)
    x1=a2-a1
    x2=a3-a1
    t1=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    x1=a1-a2
    x2=a3-a2
    t2=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    x1=a1-a3
    x2=a2-a3
    t3=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    if t1>0.0000001  && t2>0.0000001 && t3>0.0000001
        k+=1
    end
end
println(k,",",n,",",k/n)

274718,1000000,0.274718

鋭角三角形となるのは約27%でした。

単位円の内部に3点を任意に取る。

3番目は，単位円の内部に3点を取ります。これも100万回作ってその割合を調べることにしました。

circle1.jl

n=1000000
k=0
for i=1:n
    b1=rand(2)
    a1=[b1[1]*cos(2*π*b1[2]),b1[1]*sin(2*π*b1[2])]
    b1=rand(2)
    a2=[b1[1]*cos(2*π*b1[2]),b1[1]*sin(2*π*b1[2])]
    b1=rand(2)
    a3=[b1[1]*cos(2*π*b1[2]),b1[1]*sin(2*π*b1[2])]
    x1=a2-a1
    x2=a3-a1
    t1=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    x1=a1-a2
    x2=a3-a2
    t2=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    x1=a1-a3
    x2=a2-a3
    t3=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    if t1>0.00000001  && t2>0.00000001 && t3>0.00000001
        k+=1
    end
end
println(k,",",n,",",k/n)

242824,1000000,0.242824

鋭角三角形となるのは約24%でした。正方形のときより，ちょっと小さいです。

単位円の円周上に3点を取ると確率が $\frac14$ になるらしい。

4番目は，単位円の円周上に3点を取ります。これも100万回作ってその割合を調べることにしました。いろんなところで $\frac14$ であることが伝えられています。

circle2.jl


n=1000000
k=0
for i=1:n
    b1=rand(1)
    a1=[cos(2*π*b1[1]),sin(2*π*b1[1])]
    b1=rand(1)
    a2=[cos(2*π*b1[1]),sin(2*π*b1[1])]
    b1=rand(1)
    a3=[cos(2*π*b1[1]),sin(2*π*b1[1])]
    x1=a2-a1
    x2=a3-a1
    t1=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    x1=a1-a2
    x2=a3-a2
    t2=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    x1=a1-a3
    x2=a2-a3
    t3=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    if t1>0.0000001  && t2>0.0000001 && t3>0.0000001
        k+=1
    end
end
println(k,",",n,",",k/n)

250651,1000000,0.250651

100万回作って，250651個鋭角三角形ができた。約25%。ピッタリだ！！！

単位円の円周を $N$ 等分して3点を選ぶとき，鋭角三角形となる確率。

先程のtwitterの @Ihi20PvdFguJC1F さんもここから $N\to\infty$ として， $\frac14$ を得ています。

$N=2n$ のとき

円周の $2n$ 個の点から3個選ぶ。 $\frac{2n(2n-1)(2n-2)}6$ 個。

鋭角三角形となるのは $\frac{n(n-1)(n-2)}3$ 個。

よって，確率は

\frac{2n(n-1)(n-2)}{2n(2n-1)(2n-2)}=\frac{n-2}{2(2n-1)}\xrightarrow[n\to \infty]{}\frac14

$N=2n+1$ のとき

円周の $2n+1$ この点から3個選ぶ。 $\frac{(2n+1)(2n)(2n-1)}6$ 個。

鋭角三角形となるのは $\frac{n(n+1)(2n+1)}6$ 個。

よって，確率は

\frac{n(n+1)(2n+1)}{(2n+1)(2n)(2n-1)}=\frac{n+1}{2(2n-1)}\xrightarrow[n\to \infty]{}\frac14

5番目は，単位円の周上を $N$ 分割して，やっぱり100万回作ってその割合を調べることにしました。今回は理論値があるので，それと比較しました。

circle3.jl

for m=3:20
n=1000000
k=0
p=0
for i=1:n
    b1=sample(1:m,3,replace=false)
    a1=[cos(2*π*b1[1]/m),sin(2*π*b1[1]/m)]
    a2=[cos(2*π*b1[2]/m),sin(2*π*b1[2]/m)]
    a3=[cos(2*π*b1[3]/m),sin(2*π*b1[3]/m)]
    
    x1=a2-a1
    x2=a3-a1
    t1=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    
    x1=a1-a2
    x2=a3-a2
    t2=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    
    x1=a1-a3
    x2=a2-a3
    t3=(x1[1]*x2[1]+x1[2]*x2[2])/(sqrt(x1[1]^2+x1[2]^2)*sqrt(x2[1]^2+x2[2]^2))
    
    if t1>0.00000001  && t2>0.00000001 && t3>0.00000001
        k+=1
    end
end
    if m%2==0
        p=(m/2-2)/(2*(m-1))
    else
        p=((m-1)/2+1)/(2*(m-2))
    end
println("m=",m,",","100万回調べた割合は",k/n,"。理論的な確率は",p)
end

m=3,100万回調べた割合は1.0。理論的な確率は1.0
m=4,100万回調べた割合は0.0。理論的な確率は0.0
m=5,100万回調べた割合は0.500154。理論的な確率は0.5
m=6,100万回調べた割合は0.100114。理論的な確率は0.1
m=7,100万回調べた割合は0.399824。理論的な確率は0.4
m=8,100万回調べた割合は0.143472。理論的な確率は0.14285714285714285
m=9,100万回調べた割合は0.356421。理論的な確率は0.35714285714285715
m=10,100万回調べた割合は0.167043。理論的な確率は0.16666666666666666
m=11,100万回調べた割合は0.332792。理論的な確率は0.3333333333333333
m=12,100万回調べた割合は0.181575。理論的な確率は0.18181818181818182
m=13,100万回調べた割合は0.318368。理論的な確率は0.3181818181818182
m=14,100万回調べた割合は0.192279。理論的な確率は0.19230769230769232
m=15,100万回調べた割合は0.307615。理論的な確率は0.3076923076923077
m=16,100万回調べた割合は0.199948。理論的な確率は0.2
m=17,100万回調べた割合は0.299358。理論的な確率は0.3
m=18,100万回調べた割合は0.205984。理論的な確率は0.20588235294117646
m=19,100万回調べた割合は0.294719。理論的な確率は0.29411764705882354
m=20,100万回調べた割合は0.210006。理論的な確率は0.21052631578947367

バッチリな感じです！！！

まとめ

最初，内積の式がちがっていて，確率が7%くらいになっていて，「おかしいな。。。。」と1日くらい思っていました。

その後，t1>0などは，t1>0.0000001などにして，微妙なものを除くことにしました。

こうやって，実験してみると，やはり，1x1の正方形の内部の27%と単位円の内部の24%の理論値を求めたいです。（調べれば，どこかに載っていそうですが。。。。）また，どうして，単位円よりも正方形の方が鋭角三角形ができやすいのかも不明です。

わかったら，報告します！！

はじめに

確率0

1x1の正方形の内部に3点を任意に取る。

単位円の内部に3点を任意に取る。

単位円の円周上に3点を取ると確率が\frac14になるらしい。

単位円の円周をN等分して3点を選ぶとき，鋭角三角形となる確率。

まとめ

Discussion

単位円の円周上に3点を取ると確率が $\frac14$ になるらしい。

単位円の円周を $N$ 等分して3点を選ぶとき，鋭角三角形となる確率。