参照SPLISHSPLASH那个库,又学习了一下它的代码,于是发现所谓perfect sampling是这样的:
(画出来结果和mzhang 助教 所说一致)
蓝色小圈代表邻居粒子
橙色小圈是当前粒子
黑色大圈是核半径
这里取粒子半径为0.025
核半径是粒子半径的四倍
只要比较各个蓝色小圆的中心点是不是在黑色大圈里面就行了
于是这就是perfect sampling(SplishSPlasH的)
画图的代码如下
import numpy as np
import math
import matplotlib.pyplot as plt
def circle(x,y,r,color='k',count=1000):
xarr=[]
yarr=[]
for i in range(count):
j = float(i)/count * 2 * np.pi
xarr.append(x+r*np.cos(j))
yarr.append(y+r*np.sin(j))
plt.plot(xarr,yarr,c=color)
hp = 0.025
h = 4* hp
diam =2*hp
xi = [0, 0]
xj =[-h,-h]
print("h=",h,"diam=",diam)
i,j=0,0
fig1 = plt.figure(num='fig1',figsize=(6,6),dpi=100,facecolor='#FFFFFF',edgecolor='#FF0000')
plt.grid()
plt.arrow(-0.15,0,0.3,0)
plt.arrow(0,-0.15,0,0.3)
circle(0,0,h)
while(xj[0]<=h):
print("")
print("")
print("loop0 ")
while(xj[1]<=h):
j=j+1
print("loop1")
print("xi=",xi,"xj=",xj)
dist2=(xi[0]-xj[0])**2 + (xi[1]-xj[1])**2
print("dist=",math.sqrt(dist2))
#画蓝色圈
circle(xj[0],xj[1],hp,'#00BFFF')
plt.scatter(xj[0],xj[1],c='k')
plt.text(xj[0],xj[1],color='r',s=f'{j}',size=8)
plt.text(xj[0]+0.005,xj[1]+0.005,color='k',s=f'{xj[0],xj[1]}',size=5)
if dist2<= h**2:
print("yes")
xj[1] += diam
xj[0] += diam
xj[1] = -h
print("end")
#画橙色圈
circle(0,0,hp,'#FF8C00')
plt.show()