Main Content

このページの内容は最新ではありません。最新版の英語を参照するには、ここをクリックします。

シストリック アレイでの CORDIC を使用したハードウェア効率に優れた QR 分解の実装

この例では、Simulink® でハードウェア効率に優れた MATLAB® コードを使用して行列の QR 分解を計算する方法を示します。

方程式系を解く、つまり QR 分解を使用して行列方程式 AX = B の最小二乗解を計算するには、RQ'B を計算します。ここで QR = A および RX = Q'B です。R は上三角行列、Q は直交行列です。QR だけが必要な場合は、B を単位行列に設定します。

この例では、CORDIC アルゴリズムを使用してギブンス回転を適用して、行列 A から R が計算されています。C = Q'B は、同じギブンス回転を適用して行列 B から計算されています。このアルゴリズムでは反復的なシフトと加算のみを使用してこれらの計算を実行します。

モデルを開く

fxpdemo_real_4x4_systolic_array_QR_model モデルを開きます。この例では、A の行数と列数が 4 でなければなりません。行列 B には 4 つの行がなければなりませんが、列数は任意です。シストリック アレイは、同じパターンに従ってどのサイズにも拡張できます。このモデルでは、固定小数点、倍精度、および単精度のデータ型がサポートされています。この実装は、実数入力行列でのみ動作します。

open_system('fxpdemo_real_4x4_systolic_array_QR_model')

独自の入力行列 AB を入力するには、対応する Enabled Subsystem ブロックのブロック パラメーターをモデルの左側に開きます。シミュレーション後、モデルは計算された出力行列 C = Q'B および R をワークスペースに返します。"transpose(Q)*B, R 4x4 Real CORDIC Systolic Array" サブシステムのブロック パラメーターで CORDIC 反復数を指定できます。入力が固定小数点の場合、CORDIC 反復数は語長より小さくなければなりません。計算の精度は反復ごとに 1 ビットずつ、入力の語長まで向上します。

アルゴリズムでの因数分解の実行方法を確認するには、"transpose(Q)*B, R 4x4 Real CORDIC Systolic Array" サブシステムのマスク内を確認します。注釈には、下対角要素をゼロ設定するために演算されている行列の行が示されます。

CORDIC を使用してギブンス回転を実行する MATLAB Function ブロックの MATLAB コードを確認するには、引き続きブロック マスク内を確認します。

QR を使用して行列方程式 Ax = B を解く

行列方程式 AX = B を解く最初のステップは RX = Q'B を計算することです。ここで、R は上三角行列、Q は直交行列で Q*R = A です。

以下の入力は倍精度浮動小数点型です。反復数は double の仮数のビット数である 52 に設定します。

NumberOfCORDICIterations = 52;
A = 2*rand(4,4)-1;
B = 2*rand(4,4)-1;

RC = Q'B を計算するには、モデルをシミュレートします。

sim fxpdemo_real_4x4_systolic_array_QR_model
R
R = 4×4

    1.5149   -0.0519    1.7292   -0.3224
         0    0.9593   -0.0259   -0.0879
         0         0    0.2565    1.0888
         0         0         0   -0.6429

C
C = 4×4

    0.5942   -0.2382    0.0676   -0.9370
   -0.8887    0.6146   -0.5758    0.3051
    0.1725    0.7339    0.5409    0.5374
    0.8540    1.1078   -0.2183   -0.5620

R および C = Q'B で後退代入すると MATLAB と同じ結果になることを確認します。

X = R\C
X = 4×4

   -7.1245  -12.1131   -0.6637    1.4236
   -0.8779    0.7572   -0.5511    0.3545
    6.3113   10.1759    0.6673   -1.6155
   -1.3283   -1.7231    0.3396    0.8741

X_should_be = A\B
X_should_be = 4×4

   -7.1245  -12.1131   -0.6637    1.4236
   -0.8779    0.7572   -0.5511    0.3545
    6.3113   10.1759    0.6673   -1.6155
   -1.3283   -1.7231    0.3396    0.8741

組み込みの MATLAB の解と CORDIC QR の解の差のノルムは小さくなります。

norm(X - X_should_be)
ans = 3.2578e-14

QR の計算

Q および R を計算するには、B を単位行列と等価に設定します。

NumberOfCORDICIterations = 52;

A = 2*rand(4,4)-1;
B = eye(size(A,1),'like',A);

sim fxpdemo_real_4x4_systolic_array_QR_model

B が単位行列と等しい場合、Q = C' になります。

Q = C';

理論的な QR 分解は QR==A であるため、計算された QR と A の差は小さくなります。

norm(Q*R - A)
ans = 2.2861e-15

QR は一意ではない

QR 分解が一意であるのは R の行と Q の列の符号までです。R の対角要素をすべて正にすると、一意の QR 分解にすることができます。

D = diag(sign(diag(R)));
Qunique = Q*D
Qunique = 4×4

   -0.3086    0.1224   -0.1033   -0.9376
   -0.6277   -0.7636   -0.0952    0.1174
   -0.5573    0.3930    0.7146    0.1559
    0.4474   -0.4975    0.6852   -0.2877

Runique = D*R
Runique = 4×4

    1.4459   -0.8090    0.1547    0.3977
         0    1.1441    0.0809   -0.2494
         0         0    0.8193    0.1894
         0         0         0    0.4836

モデルから計算された QR と MATLAB 関数 qr を比較します。

[Q0,R0] = qr(A);
D0 = diag(sign(diag(R0)));
Q0 = Q0*D0
Q0 = 4×4

   -0.3086    0.1224   -0.1033   -0.9376
   -0.6277   -0.7636   -0.0952    0.1174
   -0.5573    0.3930    0.7146    0.1559
    0.4474   -0.4975    0.6852   -0.2877

R0 = D0*R0
R0 = 4×4

    1.4459   -0.8090    0.1547    0.3977
         0    1.1441    0.0809   -0.2494
         0         0    0.8193    0.1894
         0         0         0    0.4836

ハードウェア効率の高い実装のための固定小数点の使用

固定小数点入力データ型を使用して、ASIC や FPGA デバイス向けの効率的な HDL コードを生成します。

A および B の 3 次元配列を作成して、モデルから多くのテスト入力を実行できます。

n_test_inputs = 100;

行列 AB の乱数入力を -11 でスケーリングします。固定小数点の語長を 18 ビット、小数部の長さを 14 ビットに設定して QR 分解の成長と CORDIC アルゴリズムの中間計算を可能にします。

word_length = 18;
fraction_length = 14;

CORDIC 反復の最高精度数は語長から 1 を引いた数です。CORDIC 反復数が word_length - 1 より小さい値に設定されると、次の Ready 信号への遅延と時間刻みが短くなりますが、精度は低くなります。生成されたコードは固定小数点の語長よりも大きいシフトをサポートしていないため、CORDIC 反復数はこれより大きい値に設定できません。

NumberOfCORDICIterations = word_length - 1
NumberOfCORDICIterations = 17

乱数テスト入力は連結されて、時間 k における入力は A(:,:,k) および B(:,:,k) になります。A および B の各要素は、-1+1 の一様な確率変数になります。

A = 2*rand(4,4,n_test_inputs)-1;

Q=C' となるように単位行列として B を選択します。

B = eye(4);
B = repmat(B,1,1,n_test_inputs);

A を固定小数点にキャストし、BA と同様にキャストします。

A = fi(A,1,word_length,fraction_length);
B = cast(B,'like',A);

モデルのシミュレーション

固定小数点の入力を使用してモデルをシミュレートします。

sim fxpdemo_real_4x4_systolic_array_QR_model

誤差の計算とプロット

誤差を計算します。

norm_error = zeros(1,size(R,3));
for k = 1:size(R,3)
    Q_times_R_minus_A = double(C(:,:,k))'*double(R(:,:,k)) - double(A(:,:,k));
    norm_error(k) = norm(Q_times_R_minus_A);
end

誤差をプロットします。誤差は 10^-3 と同程度になります。

clf
plot(norm_error,'o-')
grid on
title('norm(QR - A)')

Figure contains an axes object. The axes object with title norm(QR - A) contains an object of type line.

%#ok<*NASGU,*NOPTS>

関連するトピック