如何有效地用查找替换函数？答案

【问题标题】：How do I efficiently replace a function with a lookup?如何有效地用查找替换函数？
【发布时间】：2016-02-05 01:51:13
【问题描述】：

我正在尝试提高在大型数据集上运行的代码的速度。我需要执行函数out = sinc(x)，其中x 是一个2048-by-37499 的双精度矩阵。这非常昂贵并且是我的程序的瓶颈（即使在 GPU 上计算时也是如此）。

我正在寻找任何可以提高此操作速度的解决方案。我希望这可以通过预先计算一个向量LookUp = sinc(y) 来实现，其中y 是向量y = min(min(x)):dy:max(max(x))，即一个跨越整个预期x 元素范围的向量。

如何从这个LookUp 向量有效地生成sinc(x) 的近似值？

我需要避免生成一个 3 维数组，因为这会消耗比我可用的内存更多的内存。

这是对 interp1 解决方案的测试：

a = -15;
b = 15;
rands = (b-a).*rand(1024,37499) + a;

sincx = -15:0.000005:15;
sincy = sinc(sincx);

tic
res1 = interp1(sincx,sincy,rands);
toc

tic
res2 = sinc(rands);
toc'

sincx = gpuArray(sincx);
sincy = gpuArray(sincy);
r = gpuArray(rands);

tic
r = interp1(sincx,sincy,r);
toc

r = gpuArray(rands);

tic
r = sinc(r);
toc

经过的时间是 0.426091 秒。
经过的时间是 0.472551 秒。
经过的时间是 0.004311 秒。
经过的时间是 0.130904 秒。

分别对应CPU interp1、CPU sinc、GPU interp1、GPU sinc

【问题讨论】：

我认为可以使用interp1 解决这个问题，但评估interp1 实际上比sinc 慢
您能提供一个示例输入供我们使用吗？
你的矩阵x中会有零吗？
截断泰勒近似值怎么样？可能会更快：wolframalpha.com/input/?i=sinc+taylor+approximation
@B.Thomas : interp1 代码最大逼近误差高达 28% ..

标签： matlab matrix vector

【解决方案1】：

不确定我是否完全理解您的问题。但是一旦你有了LookUp = sinc(y)，你就可以使用Matlab函数interp1

out = interp1(y,LookUp,x)

其中 x 可以是任意大小的矩阵

【讨论】：

谢谢 - 这就是我要找的！它确实会导致显着的加速（取决于查找矩阵中有多少样本点）。为了改进这一点，我将优化查找表断点的位置，以最大限度地减少特定最大错误所需的条目数。
请注意，加速越大，近似的质量就越低
@B.Thomas 我对你的实现很好奇。我对此进行了测试，发现它的速度是sinc 的两倍。在我的测试中，我使用了x = rand(2048,37499); 形式的随机数，这似乎与LookUp 中的点数无关。我确实是从工作区而不是在 GPU 上运行的。
@Matt 请查看我对上述答案的评论。对我来说，如果查找表中的样本数量很少（小于大约 6000），则 interp1() 在 CPU 上的速度会稍快一些。较大的样本数有利于 sinc() 函数。在 GPU 上，改变样本数量几乎没有影响，但 interp1() 往往比 sinc() 快 30 倍左右。

【解决方案2】：

m = min(x(:));
y = m:dy:max(x(:));
LookUp = sinc(y);

现在sinc(n) 应该等于

LookUp((n-m)/dy + 1)

假设n 是dy 的整数倍，并且位于m 和max(x(:)) 范围内。要获得LookUp 索引（即1 和numel(y) 之间的整数，我们首先移动n 但最小m，然后将其缩放dy，最后添加1，因为MATLAB 索引来自1 而不是0。

我不知道这对你的效率有什么影响，但试试看吧。

您也可以将其放入匿名函数中以提高可读性：

sinc_lookup = @(n)(LookUp((n-m)/dy + 1))

现在你可以打电话了

sinc_lookup(n)

【讨论】：

【解决方案3】：

我得出的结论是，您的代码无法显着改进。最快的查找表是基于简单的索引。对于性能测试，我们只根据随机数据进行测试：

%test data:
x=rand(2048,37499);
%relevant code:
out = sinc(x);

现在基于整数索引的查找：

a=min(x(:));
b=max(x(:));
n=1000;
x2=round((x-a)/(b-a)*(n-1)+1);
lookup=sinc(1:n);
out2=lookup(x2);

无论查找表或输入数据的大小如何，两个代码块中的最后一行所花费的时间大致相同。让sinc 评估的速度大致与索引操作一样快，我只能假设它已经使用查找表实现了。

【讨论】：

在 sinc 评估中需要时间的是调用find(x==0)。如果 OP 在他的矩阵中没有零，则重新实现的速度可以提高 3 倍
@BillBokeey：这很好，使用sin(pi*x)./(pi*x) 它比索引更快。我尝试做一些详细的基准测试，但我的 matlab 版本不断遇到该数据大小的内存不足异常。
你想要做哪些基准测试@Daniel？我认为我的电脑可以处理这个
@BillBokeey：我比较了out = sinc(x);、out2=lookup(x2); 和out3 = sin(pi*x)./(pi*x);。我正在使用一个旧的 32 位版本的 matlab，如果你有一个最近的版本，它的功能是通过timeit 运行这三行
我不太确定如何使用timeit（比如：timeit(@() sinc(x))？）

【解决方案4】：

我找到了一种更快的方法（如果您的 PC 上有 NVIDIA GPU），但是对于 x=0，这将返回 NaN，但是如果出于任何原因，您可以处理 NaN 或者您知道那么它永远不会为零：

如果您定义 r = gpuArray(rands); 并实际在 GPU 中自行评估 sinc 函数为：

tic
r=rdivide(sin(pi*r),pi*r);
toc

这通常比 GPU 中的interp1 版本的速度快 3.2 倍，而且更准确（使用上面的代码进行测试，使用不同的随机数据迭代 100 次，两种方法都具有相似的标准）。

这是因为sin 和元素除法rdivide 也是GPU 实现的（而由于某种原因sinc 不是）。见：http://uk.mathworks.com/help/distcomp/run-built-in-functions-on-a-gpu.html

【讨论】：

请注意，这仅在您拥有支持 CUDA 的显卡（即 Nvidia 显卡）时才有效。