5

I am using the matlab code from this book: http://books.google.com/books/about/Probability_Markov_chains_queues_and_sim.html?id=HdAQdzAjl60C Here is the Code:

    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);

Is there a faster code that I am not aware of? I am calling this functions millions of times, and it takes too much time.

4

7 に答える 7

9

MATLAB には組み込みの線形代数ルーチンの完全なセットがhelp slashありhelp luますhelp chol。

内部では、これらの関数は通常、最適化されたLAPACK/BLASライブラリ ルーチンを呼び出しています。これは一般に、どのプログラミング言語でも線形代数を実行する最速の方法です。MATLAB のような「遅い」言語と比較すると、m ファイルの実装より桁違いに高速であっても、予想外ではありません。

お役に立てれば。

于 2012-01-27T04:02:55.530 に答える
8

特に独自の実装を検討している場合を除き、Matlab のバックスラッシュ演算子 ( mldivide) を使用するか、係数が必要な場合はlu. mldivide は、ガウス消去法以上のことを実行できることに注意してください (たとえば、適切な場合、線形最小二乗法を実行します)。

mldivideおよびで使用されるアルゴリズムluは C および Fortran ライブラリからのものであり、Matlab での独自の実装は決して高速ではありません。ただし、独自の実装を使用することを決定し、それをより高速にしたい場合は、実装をベクトル化する方法を探すことが 1 つのオプションです (ここから始めることもできます)。

注意すべきもう1つのこと：質問の実装はピボットを行わないため、数値安定性は通常、ピボットを行う実装よりも悪くなり、一部の非特異行列では失敗することさえあります.

ガウス消去法のさまざまなバリエーションが存在しますが、それらはすべて O( n 3 ) アルゴリズムです。あるアプローチが別のアプローチよりも優れているかどうかは、特定の状況によって異なり、さらに調査する必要があります。

于 2012-01-27T04:04:58.147 に答える
0

n 行 n 列の単純なアプローチ (別名、行スワップなし) の場合:

function A = naiveGauss(A)

% find's the size

n = size(A);
n = n(1);

B = zeros(n,1);

% We have 3 steps for a 4x4 matrix so we have
% n-1 steps for an nxn matrix
for k = 1 : n-1     
    for i = k+1 : n
        % step 1: Create multiples that would make the top left 1
        % printf("multi = %d / %d\n", A(i,k), A(k,k), A(i,k)/A(k,k) )
        for j = k : n
            A(i,j) = A(i,j) - (A(i,k)/A(k,k)) *  A(k,j);
        end
        B(i) = B(i) - (A(i,k)/A(k,k))  * B(k);
    end
end
于 2012-09-27T03:07:12.150 に答える
0

matlab 関数 rref を使用できると思います。

[R,jb] = rref(A,tol)

縮小された行エシェロン形式でマトリックスを生成します。私の場合、それは最速の解決策ではありませんでした。以下のソリューションは、私の場合、約 30% 高速でした。

function C = gauss_elimination(A,B)
i = 1; % loop variable
X = [ A B ]; 
[ nX mX ] = size( X); % determining the size of matrix  

while i <= nX % start of loop 
    if X(i,i) == 0 % checking if the diagonal elements are zero or not
        disp('Diagonal element zero') % displaying the result if there exists zero 
        return
    end
    X = elimination(X,i,i); % proceeding forward if diagonal elements are non-zero
    i = i +1;
end

C = X(:,mX);

function X = elimination(X,i,j)
% Pivoting (i,j) element of matrix X and eliminating other column
% elements to zero
[ nX mX ] = size( X); 
a = X(i,j);
X(i,:) = X(i,:)/a;

for k =  1:nX % loop to find triangular form 
    if k == i
        continue
    end
    X(k,:) = X(k,:) - X(i,:)*X(k,j); 
end
于 2019-11-04T18:03:47.240 に答える