时间戳

首页 > 解决方案 > 计算 CFD 网格中点之间的距离

python - 计算 CFD 网格中点之间的距离

问题描述

我有一个大约 3e+6 点的结构化网格，并且

请查看图片： https ://i.stack.imgur.com/eUUkQ.jpg

物理域（欧几里得）中的每个点在计算域中都有一个索引（i、j 和 k）。

我需要遍历计算域索引并使用对应点进行计算。

例如，给定索引处i方向的长度将是（伪代码）： length = vec_len(point( i+1, j, k ) - point( i, j, k ))

找出给定点的长度比也很重要。例如，我将计算i方向上两个附近的长度并将它们分开。

我想出的东西花费了太多时间，并且可能没有充分利用 NumPy 所提供的全部潜力。

我制作了一个填充零的 ndarray，它将保存所有网格 XYZ 坐标。

block_data =numpy.zeros((i_dim,  j_dim, k_dim, 3), dtype='float')

数字 3 对应于 3 个元素，x、y 和 z。

因此，如果我想要i=3、j=7、k=10处的 z 值，它将是：

Z = block_data[3][7][10][2]

欧几里得空间中的一个点将是一个 (1,3) ndarray：

point = block_data[i][j][k]

我计算两点之间长度的方法是：

numpy.linalg.norm(point2 - point1)

只有长度计算大约需要 1.5 毫秒，我想计算所有点和所有方向的距离：3e+6 * 3。

我认为构造主块 ndarray (block_data) 的方法存在问题，因为它限制了我一次只能对两个点进行计算，即只有两个小的 (1,3) ndarray。如果我没记错的话，对小数组进行计算并不是那么有效。

我怎样才能解决问题并使运行时间更快？有没有针对这类问题的书推荐？谢谢 ：-）

标签： pythonnumpy

解决方案

要计算这种大小的数组中的欧几里德距离，我建议使用矢量化方法：

def euclid_dist(array, direction):
    if direction == 'i':  # make shifted views depending on the direction
        p1 = array[:-1, :, :]
        p2 = array[1:, :, :]
    elif direction == 'j':
        p1 = array[:, :-1, :]
        p2 = array[:, 1:, :]
    elif direction == 'k':
        p1 = array[:, :, :-1]
        p2 = array[:, :, 1:]
    else:
        raise ValueError('direction ' + direction + ' not known.')

    # get euclidean distance for all points in direction:
    euc_dist = (((p1 - p2)*(p1 - p2)).sum(axis=3))**0.5

    return euc_dist

使用一个小的测试数组：

arr = np.random.randint(-20, 20, 5*5*5*3).reshape(5, 5, 5, 3)
eu_i = euclid_dist(arr, 'i')
eu_j = euclid_dist(arr, 'j')
# test some values:
print(eu_i[2, 1, 2] == np.linalg.norm(arr[2, 1, 2] - arr[3, 1, 2]))
# Out 64: True
print(eu_j[1, 1, 1] == np.linalg.norm(arr[1, 1, 1] - arr[1, 2, 1]))
# Out 65: True

8e6带有点和24e6值的大数组的一些时间安排：

big_arr = np.random.rand(200, 200, 200, 3)
%timeit euclid_dist(big_arr, 'i')
# 644 ms ± 57.9 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

恕我直言，对于这样大小的数组来说，这相当快。:) 如果我正确阅读了您的时间，这比您的代码快 19000 倍。

推荐阅读