自依赖 MATLAB 循环的向量化

Question

我对 matlab 还很陌生，正在尝试理解矢量化。我通常能够找到一种方法来矢量化我的代码，但这次不同。请注意， e() 和 f() 不是函数调用，而是数组， a(i,j) 是一个矩阵。这里的问题是 e(i) 的值取决于 e(i-1) 的值。同样的事情也发生在 f(i) 上。有没有办法对这个循环进行矢量化以提高速度？

该函数采用 A 三对角对称矩阵，b 是 Ax=b 线性系统的解，并输出 sol：系统的解。代码按原样工作，我只想在可能的情况下让它更快。希望其他详细信息有所帮助。

function [sol] = MySolTridiagDirect(A,b)
   %Adds 2 slack variables
   a = [zeros(size(A,1),1),A,zeros(size(A,1),1)];
   %Adds solutions to slack variables 
   b = [zeros(1,size(b,2));b;zeros(1,size(b,2))];

   e = NaN(size(A,1)+1,1);
   f = NaN(size(A,1)+1,1);
   x = NaN(size(A,1)+1,1);

   e(1) = 0;
   f(1) = 0;
   x(1) = 0;

   for i=2:size(f,1)
       e(i) = (-a(i-1,i+1)) / (a(i-1,i-1)*e(i-1,1) + a(i-1,i));
       f(i) = (b(i) - a(i-1,i-1)*f(i-1)) / (a(i-1,i-1)*e(i-1,1) + a(i-1,i));    
   end
   %% Solver for variable 'x' (solution)
   x(end) = f(end);
   for i=size(f,1)-1:-1:2
      x(i) = e(i)*x(i+1) + f(i);
   end
   sol = x(2:end,:);
end

Answer 1

我认为这个问题不能以传统方式向量化，我没有发现任何简化代码的可能性，这样做可以将性能提高约 2 倍

a2=diag(a,2);
a0=diag(a,0);
a1=diag(a,1);

for i=2:size(f,1)
    e(i) = (-a2(i-1)) / (a0(i-1)*e(i-1) + a1(i-1));
    f(i) = (b(i) - a0(i-1)*f(i-1)) / (a0(i-1)*e(i-1) + a1(i-1));    
end

为了找到矢量化，我使用符号工具箱显示完整的公式来仅使用 e(1) 来计算 e(i)，但我看不到以矢量化方式实现这些计算的任何可能性。如果您查看公式，则涉及的嵌套括号越来越多：

n=4;
f=sym('f',[n,1]);
e=sym('e',[n,1]);
a1=sym('a1',[n,1]);
a0=sym('a0',[n,1]);
a2=sym('a2',[n,1]);
b=sym('b',[n,1]);

for i=2:size(f,1)
    e(i) = simplify((-a2(i-1)) / (a0(i-1)*e(i-1) + a1(i-1)));
    f(i) = simplify((b(i) - a0(i-1)*f(i-1)) / (a0(i-1)*e(i-1) + a1(i-1)));    
end

我知道评估此类公式的唯一可能性是递归和 for 循环。由于 for 循环比 matlab 中的递归快得多，因此可能没有更好的方法。