Matlab中的高斯消元法

Question

我使用的是本书中的matlab代码:http://books.google.com/books/about/Probability_Markov_chains_queues_and_sim.html？id = HDAQdzAjl60C 以下是代码:

    function [pi] = GE(Q)

    A = Q';
    n = size(A);
    for i=1:n-1
      for j=i+1:n
         A(j,i) = -A(j,i)/A(i,i);
      end
         for j =i+1:n
            for k=i+1:n
        A(j,k) = A(j,k)+ A(j,i) * A(i,k);
         end
      end
      end

     x(n) = 1;
      for i = n-1:-1:1
        for j= i+1:n
          x(i) = x(i) + A(i,j)*x(j);
        end
       x(i) = -x(i)/A(i,i);
      end

      pi = x/norm(x,1);

是否有一个我不知道的更快的代码？我称这个功能数百万次,而且需要花费太多时间.

Answer 1

MATLAB有一整套内置的线性代数例程 - 类型help slash,help lu或者help chol开始介绍一些在MATLAB中有效求解线性方程的常用方法.

在引擎盖下,这些函数通常调用优化LAPACK/ BLAS库例程,这通常是在任何编程语言中进行线性代数的最快方法.与像MATLAB这样的"慢"语言相比,如果它们比m文件实现快几个数量级,则不会出乎意料.

希望这可以帮助.

Answer 2

除非你特别想要实现自己的,否则你应该使用Matlab的反斜杠运算符(mldivide),或者,如果你想要这些因素,lu.请注意,mldivide可以做的不仅仅是高斯消除(例如,在适当时它会做线性最小二乘).

使用的算法通过mldivide和lu来自C和Fortran库,并在Matlab您自己的实现绝不会那么快.但是,如果您决定使用自己的实现并希望它更快,那么可以选择一种方法来实现矢量化(也许从这里开始).

还有一点需要注意:问题的实现没有进行任何转动,因此它的数值稳定性通常比转动的实现更糟糕,甚至对于某些非奇异矩阵也会失败.

存在高斯消元的不同变体,但它们都是O(n ³)算法.如果任何一种方法比另一种更好,取决于您的具体情况,您需要进行更多调查.

Answer 3

function x = naiv_gauss(A,b);
n = length(b); x = zeros(n,1);
for k=1:n-1 % forward elimination
      for i=k+1:n
           xmult = A(i,k)/A(k,k);
           for j=k+1:n
             A(i,j) = A(i,j)-xmult*A(k,j);
           end
           b(i) = b(i)-xmult*b(k);
      end
end
 % back substitution
x(n) = b(n)/A(n,n);
for i=n-1:-1:1
   sum = b(i);
   for j=i+1:n
     sum = sum-A(i,j)*x(j);
   end
   x(i) = sum/A(i,i);
end
end