thomson.py 文件源码

def ang_potential(x0):
    """
    If distance is computed along sphere rather than through 3-space.
    """
    theta = x0[0:x0.size/2]
    phi = np.pi/2-x0[x0.size/2:]

    indices = np.triu_indices(theta.size, k=1)

    theta_i = np.tile(theta, (theta.size, 1))
    theta_j = theta_i.T
    phi_i = np.tile(phi, (phi.size, 1))
    phi_j = phi_i.T
    d = _angularSeparation(theta_i[indices], phi_i[indices], theta_j[indices], phi_j[indices])
    U = np.sum(1./d)
    return U