python，向量之间的角度：需要计算效率答案

【问题标题】：python, angle between vectors: need for computational efficiencypython，向量之间的角度：需要计算效率
【发布时间】：2019-06-12 09:56:09
【问题描述】：

我有超过 500 万对两个 3D 矢量，我需要计算每对之间的角度。我试过了：

# calculate angle 
def unit_vector(vector):
    return vector / np.linalg.norm(vector)

def angle_between(v1, v2):
    v1_u = unit_vector(v1)
    v2_u = unit_vector(v2)
    return np.arccos(np.clip(np.dot(v1_u, v2_u), -1.0, 1.0))

如this post

例如：

a=[[1,1,1],[1,2,1],[6,4,5]]
b=[[1,0,0],[6,2,2],[1,9,2]]

anglevec=np.zeros(len(a))

for i in range(len(a)):
    anglevec[i]=angle_between(a[i], b[i])

print(anglevec)

但是使用循环的实现太慢了。

谁能帮忙？

【问题讨论】：

问题是：你真的需要角度吗？通过许多 3d（和 2d）计算，您通常可以跳过角度 (arccos) 并比较点积（或 products）本身。 Arccos 是这里性能最差的，点积是由简单的计算组成的。
如果您足够绝望，您可以创建一个 arccos 查找表，例如包含数千个条目。然后可以将参数缩放并截断为整数，即表中的索引。查找速度很快。你的精确度会受到限制，但你可以估计错误，看看你是否能忍受它。
当我预先计算单位向量时，我可以在不循环 a 和 b 的情况下进行点积吗？

标签： python performance loops math angle

【解决方案1】：

Numba 方法

这个任务很容易实现，例如。 https://stackoverflow.com/a/16544330/4045774 唯一的问题是 Python 循环非常慢。这可以通过使用 Numba 或 Cython 来避免。

示例

import numba as nb
import numpy as np

#You can disable parallelization with parallel=False
@nb.njit(fastmath=True,error_model="numpy",parallel=True)
def angle(v1,v2):
    #Check the dimensions, this may also have an effect on SIMD-vectorization
    assert v1.shape[1]==3
    assert v2.shape[1]==3
    res=np.empty(v1.shape[0])

    for i in nb.prange(v1.shape[0]):
        dot=0.
        a=0.
        b=0.
        for j in range(3):
            dot+=v1[i,j]*v2[i,j]
            a+=v1[i,j]**2
            b+=v2[i,j]**2
        res[i]=np.arccos(dot/(np.sqrt(a*b)))
    return res

时间

#Use numpy-arrays when working with arrays
a=np.random.rand(500_000,3)
b=np.random.rand(500_000,3)

%timeit res=your_func(a,b)
5.65 s ± 45.5 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
%timeit res=angle(a,b) #without using Numba 
3.15 s ± 13 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
%timeit res=angle(a,b) #with Numba
2.13 ms ± 441 µs per loop (mean ± std. dev. of 7 runs, 1 loop each)

【讨论】：

【解决方案2】：

这里是使用 Pool 和 Pool.map 的解决方案。

from multiprocessing import Pool, cpu_count
import numpy as np

def unit_vector(vector):
    return vector / np.linalg.norm(vector)

def angle_between(v1, v2):
    #if you want to maximize speed, avoid making variables
    #v1_u = unit_vector(v1)
    #v2_u = unit_vector(v2)
    return np.arccos(np.clip(np.dot(unit_vector(v1_u),unit_vector(v2_u)),-1.0,1.0))

def calc_angle(_list):
    #use list unpacking instead of instantiate variables
    return angle_between(*_list) 


a=[[1,1,1],[1,2,1],[6,4,5]]
b=[[1,0,0],[6,2,2],[1,9,2]]

with Pool(cpu_count()) as p: #use the context manager
    angles = p.map(calc_angle, zip(a,b))

输出：

>>> angles
[0.9553166181245092, 0.7398807743787404, 0.8775836986593762]

【讨论】：