性能 – 在矩阵的每一行上应用函数的最快方法是什么？

到目前为止我已经搜索过,我知道有几种方法( 1, 2, 3, 4)到目前为止我使用了以下代码：

Fv_calc(:,2) = arrayfun(@(n) MaxPositiveRoot2DegreePolynomial(QuadraticCoefficients(n,:)), 1:size(QuadraticCoefficients,1));

其中MaxPositiveRoot2DegreePolynomial是以下函数：

function root = MaxPositiveRoot2DegreePolynomial(c)
    d = c./c(1);
    root = eig([0 -d(3);1 -d(2)]);
    condition = root == abs(root);
    root = root(condition);
    if isempty(root)
        root = 0;
    end
    root = max(root);
end

其中,QuadraticCoefficients是62271×3矩阵,每行包含一般二次方程的系数a,b,c. ax ^ 2 bx c
关于我正在解决的问题,所有的根都是真实的,因此我使用修改的函数来查找根,以便不浪费时间检查根是否真实.
通过分析代码,我发现函数MaxPositiveRoot2DegreePolynomial的每次运行大约需要0.000047秒,对于62271运行大约需要2.914925371秒.但是代码

Fv_calc(:,2) = arrayfun(@(n) MaxPositiveRoot2DegreePolynomial(QuadraticCoefficients(n,:)), 1:size(QuadraticCoefficients,1));

需要3.198597803才能运行.所以大约0.283672431就是由arrayfun完成的.我想看看有没有办法减少这个时间？

这是一个完美的例子,它显示了你应该继续使用循环的地方.我想要的矢量化解决方案最终比你的循环方法慢.

但它可能依赖于你的Matlab版本,因为我已经将Matlab R2015b与optimzed JIT-Compiler一起使用,在旧版本中,矢量化方法可能更快.

[n,~] = size(C);
%// division
d = bsxfun(@rdivide,C,C(:,1));
%// create blocks for eigen matrix block procession
X = ( [zeros(n,1) ones(n,1) -d(:,3) -d(:,2)] ).'; %'
Y = reshape(X,2,[]);
%// use block processing to get eigenvalues of blocks
rt = blockproc(Y,[2 2],@(x) eig(x.data) ).'; %'
%// check if eigen values are real positives
condition = rt == abs(rt);
%// output filter
out = rt(condition);

基准

function [t] = bench()
    %// load your data
    DATA = load('matlab');
    C = DATA.QuadraticCoefficients;

    % functions to compare
    fcns = {
        @() thewaywewalk(C);
        @() sepideh(C);
    };

    % timeit
    t = zeros(2,1);
    for ii = 1:10;
        t = t + cellfun(@timeit, fcns);
    end
end

function out = thewaywewalk(C)  %thewaywewalk
    [n,~] = size(C);
    d = bsxfun(@rdivide,C,C(:,1));
    X = ( [zeros(n,1) ones(n,1) -d(:,3) -d(:,2)] ).'; %'
    Y = reshape(X,2,[]);
    rt = blockproc(Y,[2 2],@(x) eig(x.data) ).'; %'
    condition = rt == abs(rt);
    out = rt(condition);
end
function out = sepideh(C)  %sepideh
    [n,~] = size(C);
    out = zeros(n,1);
    for ii = 1:n
        c = C(ii,:);
        d = c./c(1);
        rt = eig([0 -d(3);1 -d(2)]);
        condition = rt == abs(rt);
        rt = rt(condition);
        if isempty(rt)
            rt = 0;
        end
        rt = max(rt);
        out(ii) = rt;
    end
end

ans =

    12.2086  %// thewaywewalk
     9.2176  %// sepideh

我仍然没有看到if条件的重点.