构建一个表示 3D 平面(存储其 Z 值)的 2D 数组,由 3 个点和所需的数组大小定义
Build a 2D array representing a 3D plane (storing its Z-values) as defined by 3 points and the desired size of the array
短
我需要找到一种优化的方法来构建表示 3D 平面的 2D numpy 数组(即存储其 Z 值)给定三个点(以数学方式定义平面)和所需的数组大小。
详情
给定三个点 p0
、p1
、p2
和一个 size
变量:
import numpy as np
import pylab as plt
p0 = [0.2, -0.4, -0.2]
p1 = [-0.6, 0.1, -0.8]
p2 = [-0.1, 0.3, -0.6]
size = (12000, 17000)
以下函数展示了如何从这 3 个点构建平面的思路:
def plane(p0, p1, p2, size):
p0 = np.array(p0)
p1 = np.array(p1)
p2 = np.array(p2)
ux, uy, uz = [p1[0] - p0[0], p1[1] - p0[1], p1[2] - p0[2]]
vx, vy, vz = [p2[0] - p0[0], p2[1] - p0[1], p2[2] - p0[2]]
uv = [uy * vz - uz * vy,
uz * vx - ux * vz,
ux * vy - uy * vx]
norm = np.array(uv)
d = -np.array(p0).dot(norm)
xx, yy = np.meshgrid(range(size[1]), range(size[0]))
z = np.array((-norm[0] * xx - norm[1] * yy - d) * 1. / norm[2])
return(z)
结果 z
如下所示:
这个实现在几秒钟内吃掉了我的记忆:
问题
您将如何优化它以尽可能以最低的内存使用量获得完全相同的结果(快速也不是可选的...)?
注意:size
在两个维度上可能大 10 倍 -> 这显然使我的 Python 在主流 2016 笔记本电脑上崩溃。
(使用最常见的科学图书馆就可以了)
我正在使用 Python 3.6.9(默认,2021 年 1 月 26 日,15:33:00)[GCC 8.4.0] Ubuntu linux 18.04.5 .
到 d
的计算为止,您的数组很小。之后的一切都是为操作的每个部分创建一个单独的 1.5GiB (12k * 17k * 8bytes/float) 临时数组。您可以完全摆脱这些数组中的大部分,并就地完成其余的工作。例如:
xx, yy = np.meshgrid(...)
您真的不需要 meshgrid:在 np.arange(size[0])[:, None]
和 np.arange(size[1])
上广播将生成您需要的所有数组。对于那部分数据,这会将您的大小从 ~3GiB 减少到 ~230KiB:
xx = np.arange(size[1], dtype=float)
yy = np.arange(size[0], dtype=float)[:, None]
表达式 z = ...
中大约有七个临时数组,每个数组的大小为 1.5GiB:
t1 = -norm[0] * xx
t2 = norm[1] * yy
t3 = t1 - t2
t4 = t3 - d
t5 = t4 * 1
t6 = t5 / norm[2]
t7 = np.array(t6)
其中一些数组可能会在计算过程中被垃圾回收,但它们都需要在某个时候进行分配。您可以使用广播和 in-place operators or the out
parameter to the appropriate ufunc
s:
摆脱所有这些
所以:
d /= norm[2]
norm /= norm[2]
xx *= -norm[0]
yy *= -norm[1]
yy -= d
z = xx + yy
即使在代码的紧凑部分,你也在做一堆不必要的操作。例如,您调用 np.array(p0)
两次。这是您的函数的更好版本:
def plane2(p0, p1, p2, size):
u = np.subtract(p1, p0)
v = np.subtract(p2, p0)
uxv = np.cross(u, v)
d = -np.dot(p0, uxv)
d /= uxv[2]
uxv /= uxv[2]
xx = np.arange(size[1], dtype=float)
xx *= -uxv[0]
yy = np.arange(size[0], dtype=float)[:, None]
yy *= -uxv[1]
yy -= d
return xx + yy
对于比建议尺寸小 10 倍的尺寸:
%timeit plane(p0, p1, p2, size)
70.1 ms ± 329 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)
%timeit plane2(p0, p1, p2, size)
8.29 ms ± 302 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
当我将大小增加到您问题中的大小时,您的实现崩溃了,而我的实现在合理的时间内完成:
%timeit plane(p0, p1, p2, size)
...
MemoryError: Unable to allocate 1.52 GiB for an array with shape (12000, 17000) and data type float64
%timeit plane2(p0, p1, p2, size)
819 ms ± 37.6 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
根据您拥有多少 RAM,以及 python 进程中的空闲内存设置,您可能无法将数组大小增加 100 倍至 152GiB,无论您执行任务的效率如何操作。
短
我需要找到一种优化的方法来构建表示 3D 平面的 2D numpy 数组(即存储其 Z 值)给定三个点(以数学方式定义平面)和所需的数组大小。
详情
给定三个点 p0
、p1
、p2
和一个 size
变量:
import numpy as np
import pylab as plt
p0 = [0.2, -0.4, -0.2]
p1 = [-0.6, 0.1, -0.8]
p2 = [-0.1, 0.3, -0.6]
size = (12000, 17000)
以下函数展示了如何从这 3 个点构建平面的思路:
def plane(p0, p1, p2, size):
p0 = np.array(p0)
p1 = np.array(p1)
p2 = np.array(p2)
ux, uy, uz = [p1[0] - p0[0], p1[1] - p0[1], p1[2] - p0[2]]
vx, vy, vz = [p2[0] - p0[0], p2[1] - p0[1], p2[2] - p0[2]]
uv = [uy * vz - uz * vy,
uz * vx - ux * vz,
ux * vy - uy * vx]
norm = np.array(uv)
d = -np.array(p0).dot(norm)
xx, yy = np.meshgrid(range(size[1]), range(size[0]))
z = np.array((-norm[0] * xx - norm[1] * yy - d) * 1. / norm[2])
return(z)
结果 z
如下所示:
这个实现在几秒钟内吃掉了我的记忆:
问题
您将如何优化它以尽可能以最低的内存使用量获得完全相同的结果(快速也不是可选的...)?
注意:size
在两个维度上可能大 10 倍 -> 这显然使我的 Python 在主流 2016 笔记本电脑上崩溃。
(使用最常见的科学图书馆就可以了)
我正在使用 Python 3.6.9(默认,2021 年 1 月 26 日,15:33:00)[GCC 8.4.0] Ubuntu linux 18.04.5 .
到 d
的计算为止,您的数组很小。之后的一切都是为操作的每个部分创建一个单独的 1.5GiB (12k * 17k * 8bytes/float) 临时数组。您可以完全摆脱这些数组中的大部分,并就地完成其余的工作。例如:
xx, yy = np.meshgrid(...)
您真的不需要 meshgrid:在 np.arange(size[0])[:, None]
和 np.arange(size[1])
上广播将生成您需要的所有数组。对于那部分数据,这会将您的大小从 ~3GiB 减少到 ~230KiB:
xx = np.arange(size[1], dtype=float)
yy = np.arange(size[0], dtype=float)[:, None]
表达式 z = ...
中大约有七个临时数组,每个数组的大小为 1.5GiB:
t1 = -norm[0] * xx
t2 = norm[1] * yy
t3 = t1 - t2
t4 = t3 - d
t5 = t4 * 1
t6 = t5 / norm[2]
t7 = np.array(t6)
其中一些数组可能会在计算过程中被垃圾回收,但它们都需要在某个时候进行分配。您可以使用广播和 in-place operators or the out
parameter to the appropriate ufunc
s:
所以:
d /= norm[2]
norm /= norm[2]
xx *= -norm[0]
yy *= -norm[1]
yy -= d
z = xx + yy
即使在代码的紧凑部分,你也在做一堆不必要的操作。例如,您调用 np.array(p0)
两次。这是您的函数的更好版本:
def plane2(p0, p1, p2, size):
u = np.subtract(p1, p0)
v = np.subtract(p2, p0)
uxv = np.cross(u, v)
d = -np.dot(p0, uxv)
d /= uxv[2]
uxv /= uxv[2]
xx = np.arange(size[1], dtype=float)
xx *= -uxv[0]
yy = np.arange(size[0], dtype=float)[:, None]
yy *= -uxv[1]
yy -= d
return xx + yy
对于比建议尺寸小 10 倍的尺寸:
%timeit plane(p0, p1, p2, size)
70.1 ms ± 329 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)
%timeit plane2(p0, p1, p2, size)
8.29 ms ± 302 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
当我将大小增加到您问题中的大小时,您的实现崩溃了,而我的实现在合理的时间内完成:
%timeit plane(p0, p1, p2, size)
...
MemoryError: Unable to allocate 1.52 GiB for an array with shape (12000, 17000) and data type float64
%timeit plane2(p0, p1, p2, size)
819 ms ± 37.6 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
根据您拥有多少 RAM,以及 python 进程中的空闲内存设置,您可能无法将数组大小增加 100 倍至 152GiB,无论您执行任务的效率如何操作。