処理の高速化

7 次查看（过去 30 天）
yuuji yamada
yuuji yamada 2019-10-30
编辑： Yoshio 2019-11-6
以前、処理の高速化ができないものかと思い、下記URLの投稿をしたものです。
https://jp.mathworks.com/matlabcentral/answers/482062-arrayfun
上記の投稿でGが単位行列と回答をしてしまいましたが私が勘違いをしておりGは
単位行列ではないことが判明しました。その場合、どのような記述になりますでしょうか。
御回答して頂いたコードにGを足して下記のように書き換えるところまではできたのですが
その先で難儀しています。どう書きかえればよいでしょうか。
%Xp = [F1(:, i), F2(:, i)];
%T(i, :, :) = ((diag(Pn(:, i)) * inv(Xp' * G * Xp) * Xp' * G) * BB)'
　　　　　　　　　　　　　　　　　　↓
%Tnew(i,:) = ((diag(Pn(:, i)) * inv(Xp' * G * Xp) * Xp' * G) * BB)'
%Tnew(i,:) = ((diag(Pn(:, i)) * ((Xp' * G * Xp)\(Xp' * G)))* BB)';
  3 个评论
显示 1更早的评论
yuuji yamada
yuuji yamada 2019-10-31
御回答ありがとうございます。
Gの構造は64ｘ64の行列になります。
Yoshio
Yoshio 2019-10-31
実対称行列であれば、ＧをLL'と分解して計算が簡略化できる可能性があるのですが、Gの例はあるでしょうか？
サンプルが無いと検証もできませんので、よろしくお願いします。

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

Yoshio
Yoshio 2019-11-4
编辑：Yoshio 2019-11-6
ベクトル化前のコードも入れていますので、参考とされてください。Gが対角行列の場合は、ベクトル化しても誤差は無視できそうですが、そうで無い場合は、誤差が大きくなる場合がありそうですので、あくまで参考にお考えください。
逆行列の公式をそのまま使うのは、数値計算上は望ましくないので、誤差評価は十分行ってください。
clear
rng
num=250000; m = 64;
%G = diag(rand(m,1));
G = rand(m);
F1=rand(m, num);
F2=rand(m, num);
Pn=rand(2, num);
BB = rand(m, 1);
T = zeros(num, 1, 2);
T0 = zeros(num,2);
Tnew1 = T0;
Tnew2 = T0;
Tnew3 = T0;
tic
for i=1:num
Xp = [F1(:, i), F2(:, i)];
%T(i, :, :) = ((diag(Pn(:, i)) * inv(Xp' * G * Xp) * Xp' * G) * BB)';
% T0(i,:) = ((diag(Pn(:, i)) * inv(Xp' * G * Xp) * Xp' * G) * BB)';
% T1(i,:) = ((diag(Pn(:, i)) * inv(Xp'*G*Xp)* Xp'*G) * BB)' ;
% T2(i,:) = ((diag(Pn(:, i)) * ((Xp'*G*Xp)\(Xp'*G)))* BB)';
% T3(i,:) = BB'*(G'*Xp*inv(Xp'*G*Xp)'*diag(Pn(:, i)));
f1 = F1(:, i);
f2 = F2(:, i);
g1 = G'*f1;
g2 = G'*f2;
A = g1'*f1;
B = g1'*f2;
C = g2'*f1;
D = g2'*f2;
Det = A*D-B*C;
a = Pn(1,i);
b = Pn(2,i);
Tnew1(i,:) = BB'*[g1 g2]*[a*[D;-B] b*[-C;A]]/Det;
Tnew2(i,:) = BB'*[a*(D*g1-B*g2) b*(-C*g1+A*g2)]/Det;
end
toc
% max(max(abs(T0-T1)))
% max(max(abs(T0-T2)))
% max(max(abs(T1-T3)))
% max(max(abs(T0-Tnew1)))
% max(max(abs(T0-Tnew2)))
max(max(abs(Tnew1-Tnew2)))
tic
G1 = G'*F1;
G2 = G'*F2;
A = sum(G1.*F1);
B = sum(G1.*F2);
C = sum(G2.*F1);
D = sum(G2.*F2);
Det = A.*D-B.*C;
Tnew3 = [Pn(1,:).*(BB'*(D.*G1-B.*G2)./Det); Pn(2,:).*(BB'*(-C.*G1+A.*G2)./Det)]';
max(max(abs(Tnew2-Tnew3)))
toc

类别

在 Help Center 和 File Exchange 中查找有关 Resizing and Reshaping Matrices 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!