Matlab中的向量化循环

Question

我正在尝试使用 Beizer 曲线通过数据点找到平滑曲线。

我有以下循环计算贝泽曲线的二项式系数：

for i = 1:n; // Here 'n' is total number of Control Points.
    a(i) = (factorial(n-1))/(factorial(i-1)*(factorial(n-i)));
end

你会如何在 Matlab 中“矢量化”这段代码？

如果我执行以下操作...我会收到错误消息。

i = 1:n;
a(i) = (factorial(n-1))/(factorial(i-1)*(factorial(n-i)));

我的另一个问题是......如果我有嵌套循环，你将如何“矢量化”它？

例如：要找到 Beizer 曲线的所有其他点...我有以下代码：

for j=2:c-1
x=0;
y=0;
for i=1:n
    kx(i)=a(i)*t(j).^(i-1)*(1-t(j)).^(n-i)*px(i);
    ky(i)=a(i)*t(j).^(i-1)*(1-t(j)).^(n-i)*py(i);
    x=x+kx(i);
    y=y+ky(i);
end
bx(j)=x;
by(j)=y;
end

Answer 1

对于您的第一个示例，有两个问题。正如@Cris Luengo 所指出的，第一个是您需要使用元素除法（以及元素乘法）。另一个半问题是您不需要使用语法 a(i)，它将尝试将整个输出数组分配到 a 的子集（使用 i 作为索引）。在这种情况下会很好，但它是不必要的，在其他一些情况下可能不会做你想要/期望的；这是一个你需要小心的习惯。在我有限的测试中，此公式产生的输出与您的原始形式相同：

i = 1:n;
a2 = (factorial(n-1))./(factorial(i-1).*(factorial(n-i)));

我应该注意，有一个内置函数 nchoosek() 基本上可以执行此计算，它支持向量输入，但仅适用于“n”项（在您的情况下不会有所不同）。因此，如果您希望将代码矢量化，它可能对您没有用处。虽然我认为更快的运行时间是您的目标，但可能值得测试一个调用 nchoosek() 的非矢量化版本，以防万一那里有一些优化巫术在起作用。您可以使用各种工具测量运行时间，例如 tic()/toc() 或分析器。

对于您的第二个示例，根据您的问题的具体情况，它可能容易或困难。您的示例依赖于足够多的预定义变量（如 c、t 等），因此对我来说测试起来并不容易。但总的来说，策略是将一个自变量定义为行向量，另一个定义为列向量。例如，假设 i 的大小为 n x 1，j 的大小为 1 x m。仅关于 i 变化的变量的大小也是 n x 1；那些仅关于 j 变化的将是 1 x m；并且两者都不同的那些将是n x m。然后，您只需要使用可以处理您提供给它们的数组的函数和运算符，并且您可能需要不使用上面矢量化示例中所需的元素运算符。这是一个有点琐碎的例子：

n = 6;
m = 4;

i = (1:n)'; % size 6 x 1
j = 1:m; % size 1 x 4

q = i * j; % size 6 x 4

所有基本的数学运算（乘法、加法等）都支持向量（通常也支持 n 维数组），并且它们会用它们做“正确的”线性代数事情。此外，根据您的 Matlab 版本，一些基本操作还支持implicit array expansion（例如，您可以向向量添加一个常数，它只会“按您的意思行事”，即使这是不正确的线性代数)

对于更复杂的情况，这可能有点难以理解。创建临时中间变量会有所帮助，这样您的代码一次做的事情就会更少。

最后，作为临别之景，我应该指出，Mathworks 历来不鼓励使用 i 和 j 作为变量名（这可能在较新的版本中有所改变）。原因是Matlab使用i和j来表示虚数。因此，如果您在代码中使用它们，那么 Matlab 必须弄清楚您的意思。这会使您的代码变慢，或者在最坏的情况下会导致错误的行为。一种常见的做法是改用 ii 和 jj，或者完全选择不同的迭代器变量。