3

我一直在努力加快以下功能，但没有结果：

function beta = beta_c(k,c,gamma)
beta = zeros(size(k));
E = @(x) (1.453*x.^4)./((1 + x.^2).^(17/6));
for ii = 1:size(k,1)
    for jj = 1:size(k,2)
        E_int = integral(E,k(ii,jj),10000);
        beta(ii,jj) = c*gamma/(k(ii,jj)*sqrt(E_int));
    end
end
end

到目前为止，我是这样解决的：

function beta = beta_calc(k,c,gamma)
k_1d = reshape(k,[1,numel(k)]);
E_1d =@(k) 1.453.*k.^4./((1 + k.^2).^(17/6));
E_int = zeros(1,numel(k_1d));
parfor ii = 1:numel(k_1d)
E_int(ii) = quad(E_1d,k_1d(ii),10000);
end
beta_1d = c*gamma./(k_1d.*sqrt(E_int));
beta = reshape(beta_1d,[size(k,1),size(k,2)]);
end

在我看来，它并没有真正提高性能。你怎么看待这件事？

你介意发光吗？

我提前谢谢你。

编辑

我将介绍一些涉及我的问题的理论背景。一般情况下，beta 计算如下

在此处输入图像描述

因此，在一维 k 数组的简化情况下，E_int 可以计算为

E = 1.453.*k.^4./((1 + k.^2).^(17/6));
E_int = 1.5 - cumtrapz(k,E);

或者，或者，作为

E_int(1) = 1.5;
for jj = 2:numel(k)
E =@(k) 1.453.*k.^4./((1 + k.^2).^(17/6));
E_int(jj) = E_int(jj - 1) - integral(E,k(jj-1),k(jj));
end

尽管如此，k目前是一个矩阵k(size1,size2)。

4

2 回答 2

2

我喜欢这个问题。

问题：该函数integral将积分限制为仅标量。因此，很难对 的计算进行矢量化E_int。

一条线索：似乎有很多冗余将相同的功能一遍又一遍地集成k(ii,jj)到无穷大......

建议的解决方案：如何将值k从最小到最大排序并E_sort_int(si) = integral( E, sortedK(si), sortedK(si+1) );与sortedK( numel(k) + 1 ) = 10000;. 然后是 的全部值E_int = cumsum( E_sort_int );（您只需要“撤消”排序并将其重新整形为 的大小k）。

于 2013-01-20T19:47:08.723 回答
2

这是另一种方法，并行化，因为它很容易使用spmdor parfor。而不是integral考虑quad，请参阅此链接以获取示例...

于 2013-01-21T07:23:30.910 回答