您可以从一些简单的相对有序的位置分布开始,然后通过应用动态系统方法/梯度体面类型迭代,您可以让位置收敛到更结构化的模式。我用python写了这样的实现,它是向量化的形式,但我还添加了一个带有for循环的等效函数,以说明函数的结构。最终的有序模式的灵感来自于稳定平衡位置,如果它们被弹簧固定,一堆相同半径 r 的圆盘会形成,每两个圆盘一个。为了简化计算,我将弹簧张力平方,从而避免平方根,因此与典型的物理模型不完全一样,但接近它。
import numpy as np
import matplotlib.pyplot as plt
def Grad(pos, r):
Dq = - pos[:, :, np.newaxis] + pos[:, np.newaxis, :]
D = Dq[0,:,:]*Dq[0,:,:] + Dq[1,:,:]*Dq[1,:,:] + np.identity(Dq.shape[1])
Dq = (1 - r**2 / D) * Dq
return - Dq.sum(axis=2)
def Grad_flow(q_, r, step):
Q = q_
n_iter = 0
while True:
n_iter = n_iter + 1 # you can count the number of iterations needed to reach the equilibrium
Q_prev = Q
Q = Q - step * Grad(Q, r)
if np.sum(np.abs((Q.T).dot(Q) - (Q_prev.T).dot(Q_prev))) < 1e-5:
return Q
'''
Test:
'''
p = np.array([[-3, 3], [-1, 3], [1,3], [3,3],
[-3, 1], [-1, 1], [1,1], [3,1],
[-3,-1], [-1,-1], [1,-1], [3,-1],
[-3,-3], [-1, -3], [1, -3], [3,-3],
[-2, 1], [-1,2],[2,-2], [-2,-2],
[2,2], [2,0]]).T
r = 0.5
step = 0.01
q = Grad_flow(p, r, step)
'''
Plot:
'''
fig, axs = plt.subplots(1,1)
axs.set_aspect('equal')
axs.plot(q[0,:], q[1,:], 'ro')
axs.plot(p[0,:], p[1,:], 'bo')
plt.grid()
plt.show()
你从蓝色位置开始,然后让它们收敛到红色位置:
这是 Grad 函数的循环版本:
def Grad(pos, r):
grad = np.zeros(pos.shape, dtype=float)
for i in range(pos.shape[1]):
for j in range(pos.shape[1]):
if not i==j:
d_pos_0 = pos[0, i] - pos[0, j]
d_pos_1 = pos[1, i] - pos[1, j]
m = d_pos_0*d_pos_0 + d_pos_1*d_pos_1
m = 1 - r*r / m
grad[0, i] = grad[0, i] + m * d_pos_0
grad[1, i] = grad[1, i] + m * d_pos_1
return grad
当然,所有这些都有点启发式,我不能保证完全通用,所以你必须玩并选择参数r,这是位置之间的一半距离,迭代步长step,初始位置p等等。