分段多项式计算的并行化

Question

我正在尝试评估从三次样条获得的大型分段多项式中的点。这需要很长时间才能完成，我想加快速度。

因此，我想用并行过程而不是按顺序计算分段多项式上的一个点。

代码：

z = zeros(1e6, 1) ;    % preallocate some memory for speed
Y = rand(11220,161) ;  %some data, rand for generating a working example
X = 0 : 0.0125 : 2 ;   % vector of data sites
pp = spline(X, Y) ;    % get the piecewise polynomial form of the cubic spline.

生成的结构很大。

for t = 1 : 1e6  % big number
    hcurrent = ppval(pp,t); %evaluate the piecewise polynomial at t
    z(t) = sum(x(t:t+M-1).*hcurrent,1) ; % do some operation of the interpolated value. Most likely not relevant to this question.
end

不幸的是，矩阵形式并使用：

hcurrent = flipud(ppval(pp, 1: 1e6 ))

需要太多内存来处理，所以无法完成。有没有办法可以批处理这段代码以加快速度？

Answer 1

对并行循环使用parfor命令。参见 here，也将 z 向量预先计算为 z(j) = x(j:j+M-1) 和 parfor 中的 hcurrent 以加快速度。

Answer 2

样条参数估计可以写成矩阵形式。

一旦你以矩阵形式编写并解决它，你就可以使用模型矩阵来评估所有数据点上的样条，使用矩阵乘法，这可能是 MATLAB 中最优化的操作。

Answer 3

对于标量第二个参数，如您的示例所示，您要处理两个问题。首先，有大量的函数调用开销和冗余计算（例如，unmkpp(pp) 被称为每个循环迭代）。其次，ppval 被编写为通用的，因此它没有完全矢量化，并且做了很多在您的情况下不必要的事情。

下面是矢量化代码代码，它利用了您问题的某些结构（例如，t 是一个大于 0 的整数），避免了函数调用开销，将一些计算移到您的主要 for 循环（以额外内存为代价），并摆脱 ppval:

内部的 for 循环

n = 1e6;
z = zeros(n,1);
X = 0:0.0125:2;
Y = rand(11220,numel(X));
pp = spline(X,Y);

[b,c,l,k,dd] = unmkpp(pp);
T = 1:n;
idx = discretize(T,[-Inf b(2:l) Inf]); % Or: [~,idx] = histc(T,[-Inf b(2:l) Inf]);
x = bsxfun(@power,T-b(idx),(k-1:-1:0).').';

idx = dd*idx;
d = 1-dd:0;
for t = T
    hcurrent = sum(bsxfun(@times,c(idx(t)+d,:),x(t,:)),2);
    z(t) = ...;
end

生成的代码占用了 n=1e6 示例的约 34% 的时间。请注意，由于矢量化，计算以不同的顺序执行。由于浮点数学的性质，这将导致 ppval 和我的优化版本之间的输出略有不同。任何差异都应该是几次 eps(hcurrent) 的数量级。您仍然可以尝试使用 parfor 来进一步加快计算速度（已经有四个运行工作人员，我的系统只用了您代码原始时间的 12%）。

我认为以上是概念验证。如果您的示例与您的实际代码和数据不符，我可能过度优化了上面的代码。在这种情况下，我建议创建您自己的优化版本。您可以通过在命令 Window 中键入 edit ppval 来查看 ppval 的代码。您可以通过查看问题的结构以及您在 z 向量中最终想要的内容来实施进一步的优化。

在内部，ppval 仍然使用 histc, which has been deprecated. My code above uses discretize to perform the same task, as suggested by the documentation。

分段多项式计算的并行化

Parallelization of Piecewise Polynomial Evaluation

multithreading

matlab

parallel-processing

vectorization

spline