3次元極座標系を用いて球面上に一様な乱数を生成しようとする。
x=r\sin{\theta}\cos{\phi}\\
y=r\sin{\theta}\sin{\phi}\\
z=r\cos{\theta}\\
\\
(0\leq{\theta} \leq {\pi}, 0\leq{\phi} \leq 2{\pi})
#実際のプログラム
今回はAnacondaのJupyter Notebookを用いてプログラムを実行した。
###全体
全体としては下記のように記述した。
import numpy as np
import random
import math
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
fig = plt.figure(figsize = (5,5))
ax = fig.add_subplot(111, projection = "3d")
test=1000
r=1.0
for i in range(test):
a = random.uniform(0, np.pi)
b = random.uniform(0, 2*np.pi)
x = r*np.sin(a)*np.cos(b)
y = r*np.sin(a)*np.sin(b)
z = r*np.cos(a)
X = [x]
Y = [y]
Z = [z]
ax.plot(X,Y,Z,marker="o")
ax.tick_params(labelbottom = False,
labelleft = False,
labelright = False,
labeltop = False)
plt.show()
###方針
方針としてはx,y,zについてr=1とし、θとΦをramdom関数を用いて擬似乱数として1000回設定し、点(x,y,z)がどのように球面状に分布するかを考える。
(果たして球面上に一様に分布するかどうかを確かめる。)
####必要なライブラリをインポート
numpy,matplotlib.pyplot,randomをインポートする。
import numpy as np
import random
import math
import matplotlib.pyplot as plt
####綺麗な3次元図形を描くための準備をする
サイズが対称な3次元図形となるために必要なコードを入れる。
from mpl_toolkits.mplot3d import Axes3D
fig = plt.figure(figsize = (5,5))
ax = fig.add_subplot(111, projection = "3d")
####for文の中で(x,y,z)を作成しプロット
r=1としtest=1000回の繰り返しを実行する。
a=θ, b=Φとしてそれぞれ0~π, 0~2πの範囲の中でランダムな数を生成する。
x,y,zにθ,Φを代入する。
このまま散布図としてプロットすることができなかったため、x,y,zの値をX,Y,Zのリストの中に入れてからプロットする。
test = 1000
r = 1.0
for i in range(test):
a = random.uniform(0, np.pi)
b = random.uniform(0, 2*np.pi)
x = r*np.sin(a)*np.cos(b)
y = r*np.sin(a)*np.sin(b)
z = r*np.cos(a)
X = [x]
Y = [y]
Z = [z]
ax.plot(X,Y,Z,marker = "o")
plt.show()
####補足
無駄な軸の目盛りを削除した。
ax.tick_params(labelbottom = False,
labelleft = False,
labelright = False,
labeltop = False)
###結果
極付近に点が集中してしまい、一様に分布しない。
つまり、極座標系を用いて球面上に一様乱数を生成することはできない。
#一様分布にならなかった理由
3次元の場合、球面上の微小面積は
dS=r^2\sin{\theta} d{\theta}d{\phi}\\\
すなわち、曲付近(θが小さい)ほど微小面積が小さく、小さな範囲に乱数が集中してしまうためである。
次回は一様乱数を球面上に作るための方法を考える。