-5

以下は、私が matlab で書いているコードの一部です。ここでは、次の単純な数学演算 [A][X]=[B] を実行します。[X] は不明です。私の場合、kの長さは約1600000です。したがって、配列の各要素のg1、g2、およびg3の値を取得するだけです。私は次のことを試しました

k31 = omega3./(d)      
k32 = omega3_2./(d)    
A   = [2,1,5;-2,-1,-5];    
X   = [g1;g2;g3];

for ii = 1:length(k31)    
    B = [k31(ii); k32(ii)];    
    X = pinv(A).*B;    
end

display(g1,g2,g3)

私は疑似逆を使用しているので、基本的に各 X の解を得ることができ、そこで編集を行いました....そして x は不明です。数学的には実行できますが、コーディングすることはできません。

また、次のように g1 g2 g3 の値を x と y でプロットするにはどうすればよいですか。scatter(x(1:end-1), y(1:end-1), 5, g2) および scatter(x(1:end-1), y(1:end-1), 5, g3)

4

2 に答える 2

1

ここでいくつかの仮定をしなければならないので、我慢してください。

私はあなたがこれをしたいと思う:

k31 = omega3./(d)      
k32 = omega3_2./(d)    
A   = [2,1,5;-2,-1,-5];    

X = cell(length(k31),1);
for ii = 1:length(k31)            
    X{ii} = A\[k31(ii); k32(ii)];    
end

invまたはの代わりにバックスラッシュ演算子を使用しますpinv。入力help slashして詳細情報を取得します。

inv通常、バックスラッシュ演算子はorよりもはるかに高速で正確ですpinv。また、はるかに柔軟です-あなたのケースは未決定です(明示的に解くことができない方程式が1つ不足しています)。その場合、バックスラッシュ演算子は最小二乗解を見つけます。

すべての結果をcell-array に保存する方法に注意してくださいX。これは、n番目のソリューションに次の方法でアクセスできることを意味します。

X{n}    % == [g1(n) g2(n) g3(n)]
于 2012-11-08T10:21:56.747 に答える
0

私の意見では、特異値分解から疑似逆数を作成する方がよいでしょう。Moore-Penrose 疑似逆関数は機能しますが、時々奇妙な結果になると思います。特に例の行列のランクが完全ではないため、最小二乗法は不安定になる可能性があります。

また; 反復ごとに逆数を計算しないでください!

以下は、A のランク不足が処理される例です。

A   = [2,1,5;-2,-1,-5];

% pseudo inverse of A
[m,n]   = size(A);
% get SVD
[U,S,V] = svd(A);
% check rank
r       = rank(S);
SR      = S(1:r,1:r);
% make complete if rank is not full
SRc     = [SR^-1 zeros(r,m-r);zeros(n-r,r) zeros(n-r,m-r)];
% create inverse
A_inv   = V*SRc*U.';

X=[];    
for i = 1:1600
    % this actually takes most of the time
    k31 = rand(1000, 1);
    k32 = rand(1000, 1);
    B   = [k31';k32'];
    % X is overwritten on every loop...
    X   = [X A_inv*B];
end

N_xy = 1000;
x    = rand(N_xy,1);
y    = rand(N_xy,1);
g1   = X(1,1:N_xy);

figure(1), clf
scatter(x,y,5,g1)

それはあなたが望むものではないので、1600000ポイントすべてをプロットしませんでした

于 2012-11-08T11:46:02.190 に答える