2

对于坐标为rho、phi、z的极点，我对m = 1: M和n = 1: N进行了双重求和：

m 和 n 的双重求和

我已经写了它的矢量化符号：

N = 10;
M = 10;
n = 1:N;
m = 1:M;

rho = 1;
phi = 1;
z = 1;

summ =  cos (n*z)  * besselj(m'-1, n*rho) * cos(m*phi)';

现在我需要重写这个函数来接受坐标rho、phi、z的向量（列） 。我尝试了arrayfun、cellfun、简单的for循环——它们对我来说太慢了。我知道“MATLAB 数组操作技巧和窍门”，但作为 MATLAB 初学者，我无法理解 repmat 和其他函数。

有人可以建议矢量化解决方案吗？

4

2 回答 2

1

我认为您的代码已经很好地矢量化了（对于n和m）。rho如果您希望此函数也接受/ phi/值的数组z，我建议您只需在 for 循环中处理值，因为我怀疑任何进一步的矢量化都会带来显着的改进（加上代码将更难阅读）。

话虽如此，在下面的代码中，我尝试{row N} * { matrix N*M } * {col M} = {scalar}通过对 BESSELJ 和 COS 函数的一次调用（我将每个行/矩阵/列放在第三维）。它们的乘法仍然在循环中完成（准确地说是 ARRAYFUN）：

%# parameters
N = 10; M = 10;
n = 1:N; m = 1:M;

num = 50;
rho = 1:num; phi = 1:num; z = 1:num;

%# straightforward FOR-loop
tic
result1 = zeros(1,num);
for i=1:num
    result1(i) = cos(n*z(i)) * besselj(m'-1, n*rho(i)) * cos(m*phi(i))';
end
toc

%# vectorized computation of the components
tic
a = cos( bsxfun(@times, n, permute(z(:),[3 2 1])) );
b = besselj(m'-1, reshape(bsxfun(@times,n,rho(:))',[],1)');             %'
b = permute(reshape(b',[length(m) length(n) length(rho)]), [2 1 3]);    %'
c = cos( bsxfun(@times, m, permute(phi(:),[3 2 1])) );
result2 = arrayfun(@(i) a(:,:,i)*b(:,:,i)*c(:,:,i)', 1:num);            %'
toc

%# make sure the two results are the same
assert( isequal(result1,result2) )

我使用TIMEIT函数做了另一个基准测试（提供更公平的计时）。结果与前面一致：

0.0062407    # elapsed time (seconds) for the my solution
0.015677     # elapsed time (seconds) for the FOR-loop solution

请注意，随着您增加输入向量的大小，这两种方法将开始具有相似的时序（在某些情况下，FOR 循环甚至会获胜）

于 2011-09-10T02:28:59.393 回答
0

您需要创建两个矩阵，例如m_，n_通过选择每个矩阵的元素，您可以获得和i,j所需的索引。mn

大多数 MATLAB 函数都接受矩阵和向量，并逐个元素地计算它们的结果。因此，要产生一个双重和，您可以并行计算和的所有元素 f(m_, n_)并将它们相加。

在您的情况下（请注意，.*运算符执行矩阵的元素乘法）

N = 10;
M = 10;
n = 1:N;
m = 1:M;

rho = 1;
phi = 1;
z = 1;

% N rows x M columns for each matrix
% n_ - all columns are identical
% m_ - all rows are identical
n_ = repmat(n', 1,  M);
m_ = repmat(m , N,  1);

element_nm =  cos (n_*z) .* besselj(m_-1, n_*rho) .* cos(m_*phi);
sum_all = sum( element_nm(:) );
于 2011-09-09T18:26:31.577 回答