MATLAB Answers

0

arrayfunの使い方について

Asked by yuuji yamada on 25 Sep 2019
Latest activity Commented on by yuuji yamada on 14 Oct 2019
現在、Matlab2018bを使用し、アプリを作成しているのですが、処理を高速化させたい部分があります。
それはforループを使用して約25万回ループしている箇所です。
arrayfun関数を使用するとforループを使用しないで書くことができるようなのですが下記のような処理の場合でも
arrayfun関数でforループを使用せず書くことができるものでしょうか。
改善したい箇所は大体、下記のような感じになっており、最終的に25万行のデータSを作成したいと思っています。
S=zeros(250000, 1, 2);
for i=1:250000
Xp = [F1(:, i), F2(:, i)];
S(i, :, :) = Fn(:, i) * Xp' * Xp * Xp' * C;
end

  4 Comments

Show 1 older comment
yuuji yamada 2019 年 10 月 1 日
下記にダミーデータで動作するコードを書いて見ました。よろしくお願いいたします。
num=25;
G=eye(64);
F1=rand(64, num);
F2=rand(64, num);
Pn=rand(2, num);
B = rand(64, 1);
T = zeros(num, 1, 2);
for i=1:num
Xp = [F1(:, i), F2(:, i)];
T(i, :, :) = ((diag(Pn(:, i)) * inv(Xp' * G * Xp) * Xp' * G) * B)';
end
Yoshio
2019 年 10 月 3 日
ダミーデータ作成ありがとうございます。
細かい話で恐縮ですが、Gが単位行列となっていますが、Xp' * G * Xpを計算しているのは、Gを単位行列以外の場合にも対応させたいということでしょうか?
yuuji yamada 2019 年 10 月 3 日
御回答ありがとうございます。
当方、該当処理部分の数式の意味を理解していなくて恐縮なのですが
Gの内容は固定で最初にeye(64)で定義した状態から変更されることはないです。
これで回答になっていますでしょうか。

Sign in to comment.

1 Answer

Yoshio
回答者: Yoshio
2019 年 10 月 8 日
編集済み: Yoshio
2019 年 10 月 8 日
 採用された回答

ダミーのコードを参照して、以下のようにしてみました。Tは3次元配列ですが、計算結果だけからだと2次元でもよさそうなので、Tnewとしてこれを求めています。
行列を求めるところを陽に公式を適用してベクトル化していますが、今の所arrayfunが使えるるコードが浮かびません。
clear
rng
num=250000; m = 64;
%G=eye(m);
F1=rand(m, num);
F2=rand(m, num);
Pn=rand(2, num);
BB = rand(m, 1);
T = zeros(num, 1, 2);
Tnew = zeros(num,2);
tic
for i=1:num
%Xp = [F1(:, i), F2(:, i)];
%T(i, :, :) = ((diag(Pn(:, i)) * inv(Xp' * G * Xp) * Xp' * G) * BB)'
%Tnew(i,:) = ((diag(Pn(:, i)) * inv(Xp' * Xp) * Xp') * BB)'
%Tnew(i,:) = ((diag(Pn(:, i)) * ((Xp'*Xp)\Xp'))* BB)';
%Tnew(i,:) = BB*Xp*((diag(Pn(:, i)) * inv(Xp' * Xp))'
f1 = F1(:, i);
f2 = F2(:, i);
A = f1'*f1;
B = f1'*f2;
C = f1'*f2;
D = f2'*f2;
E = [D -B;-C A]/(A*D-B*C);
a = Pn(1,i);
b = Pn(2,i);
Tnew(i,:) = BB'*[f1 f2]*([a*E(1,:);b*E(2,:)])';
end
toc

  2 Comments

Yoshio
2019 年 10 月 14 日
上記プログラムでループの中を見直したところ、Eの逆行列をベクトルで明示的に表すことで、ループ無しで書けました。25万行計算、2.3 GHz Intel Core i5のMac Miniでと0.4秒程度です。
tic
A = sum(F1.*F1);
B = sum(F1.*F2);
C = B;
D = sum(F2.*F2);
% E = [D -B;-C A]/(A*D-B*C);
Det = A.*D-B.*C;
Tnew3 = [Pn(1,:).*(BB'*(D.*F1-C.*F2)./Det); Pn(2,:).*(BB'*(-B.*F1+A.*F2)./Det)]';
max(max(abs(Tnew-Tnew3)))
toc
yuuji yamada 2019 年 10 月 14 日
御回答ありがとうございます。
さっそく試してみます。ありがとうございました。

Sign in to comment.