各対角番号に対する対角和からなる配列

9 views (last 30 days)
Yamada Taro
Yamada Taro on 5 Sep 2022
Commented: Yamada Taro on 3 Oct 2022
写真のような正方行列Aを持っている際に、自身より右下の第k対角成分の和をとった行列Bを作成したいのですが、何か簡単に書く方法はないでしょうか。
現在は縮小行列を作り、各縮小行列の第j番目の対角成分を取得する形で計算していますが、かなり遅いし効率も悪いなと感じています。
A = reshape(1:9,3,3);
B = zeros(3,3);
for i = 1:3
A_partial = A(i:end,i:end);
for j = 1:3-i+1
B(i,i+j-1) = sum(diag(A_partial,j-1));
B(i+j-1,i) = sum(diag(A_partial,-1*(j-1)));
end
end
A
B
  2 Comments
Yamada Taro
Yamada Taro on 3 Oct 2022
@Atsushi Uenoさん、有難うございます!!こんなサイトがあるんですね、見てみます。

Sign in to comment.

Accepted Answer

NO MIYA
NO MIYA on 5 Sep 2022
Edited: NO MIYA on 5 Sep 2022
もっと簡単な方法があると思いますが、とりあえずこんなのでどうでしょうか。
A⇒AA、B⇒BBと考えてください。
一応n×n(3<n)の正方行列に対応できるようにもしておきました。
padarrayを使うにはImage Processing Toolboxが必要なのですが、やっていることは配列の0埋めなので自作関数でいくらでもやりようがあると思います(結果、計算量が増えるかもしれませんが)
for文を回す回数が投稿者様のが二重ループで3(n)回以上なのに対し、下記プログラムは2(n-1)回で面倒な対角要素を省いたので最低限の改良はできたかと思います。
AA = reshape(1:9,3,3); % 基の正方行列
BB = AA; % 出力したい配列
n = size(AA,1); % 正方行列の1辺の長さ
%%% 処理
for i = 2:1:n
AA_plus = padarray(AA(i:end,i:end),[i-1,i-1],0,'post');
BB = BB + AA_plus;
end
BB
BB = 3×3
15 12 7 8 14 8 3 6 9
  6 Comments
Yamada Taro
Yamada Taro on 22 Sep 2022
ご返信ありがとうございます、こちら返事が遅くなってしまい大変申し訳ございません。
fortranは触ったことはないのですが、matlabであればできる限り既にかなり高速化がされている行列演算にできる限り頼る形にした方が良いとのこと、勉強になりました。
また速度重視ならコンパイルすべき、とのご指摘もありがとうございます。こちら試したことがないので、高速化が必要になった際に挑戦してみようと思います。ありがとうございました。

Sign in to comment.

More Answers (1)

Hernia Baby
Hernia Baby on 6 Sep 2022
Edited: Hernia Baby on 6 Sep 2022
せっかくなので@Yamada Taroさんのイメージどおりに動くよう、対角成分をうまく使う関数を作成しました
clc,clear;
A = reshape(1:9,3,3);
B = zeros(size(A));
キモとなる関数です
  1. diagを使って対角成分のみをとる
  2. flipして反転
  3. cumsumで累積和
  4. flipして順序を戻す
  5. diagで対角行列をつくる
f = @(X,x) diag(flip(cumsum(flip(diag(X,x)))),x);
上記の関数を - (行番号 - 1) ~ 列番号 - 1 分行います
-1としているのは角っこは計算する必要がないからです。
時間も計測してみましょうかね
tic
for ii = -height(A)+1:width(A)-1
B = B + f(A,ii);
end
toc
Elapsed time is 0.007458 seconds.
確認してみましょう
B
B = 3×3
15 12 7 8 14 8 3 6 9
  1 Comment
Yamada Taro
Yamada Taro on 22 Sep 2022
有難うございます。おかげさまで無事実装できました。

Sign in to comment.

Tags

Community Treasure Hunt

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

Start Hunting!